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Chapter 1 



Introduction 



The development of quantum mechanics was connected with one of the greatest conceptual leaps 
in theoretical physics. In order to describe properly microscopic phenomena, it was necessary to 
abandon the classical notion of a physical property, and to resort to a completely new formalism 
representing microobjects. This formalism lies at the heart of quantum mechanics. It makes a 
clear distinction between the state of a physical system, and the performed observation. The 
state is characterized by a wave function, whose time evolution is governed by an appropriate 
equation of motion (e.g. the Schrodinger equation for a nonrelativistic particle). If we want to 
relate the state to quantities observed in an experiment, we need to use the second element of the 
quantum mechanical formalism, i.e. the representation of the measuring apparatus in terms of 
operators acting on the wave function. According to the Born interpretation, quantum mechanics 
provides probabilistic predictions concerning the behaviour of a quantum system [|l]]. Quantum 
description of an experiment specifies only the chance that the measurement will yield a given 
outcome. In order to verify such predictions, we need to repeat the measurement many times on 
identically prepared systems, and then to compare the histogram of experimental outcomes with 
the probability distribution calculated from the theory. The result of a single measurement cannot 
be described in a deterministic way. This randomness seems to be a very fundamental feature 
of the microworld. So far, all attempts to introduce deterministic description of microscopic 
phenomena have failed, and experiments have ruled out whole classes of theories alternative to 
quantum mechanics. 

Pictorially speaking, the click on a measuring apparatus is only a faint shadow of the quantum 
state, which stays hidden behind the scene, though being the main actor. However, one might 
ask if it is possible to reveal experimentally the complete information on the state of a quantum 
system. This cannot be the case if we are given only a single copy of the system. Any measure- 
ment consists in a certain kind of interaction between the system and the detecting apparatus. 
This interaction maps some properties of the measured system onto the state of the apparatus, 
and makes them accessible to our cognition. After the measurement, the state of the system is 
perturbed by the interaction with the detector, and it is no longer described by the original wave 
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function. For example, observation of the position of a particle by scattering photons inevitably 
modifies its momentum. The deleterious character of quantum measurement was realized very 
early in the development of quantum mechanics, and it is closely related to the Heisenberg un- 
certainty principle [^. Thus, our system can be detected only once, and we cannot gain more 
information by repeating the measurement. One might try to circumvent this difficulty by de- 
signing a single measurement that would yield complete information on the quantum state. Such 
a measurement would map each possible state of the system onto a different, fully distinguish- 
able state of the apparatus. This means that all the final states of the apparatus would have to be 
mutually orthogonal, even these corresponding to nonorthogonal initial states of the measured 
system. Of course, this would violate unitarity of quantum evolution, and such a measurement 
does not exist. One could also consider cloning of quantum states, i.e. an operation that would 
generate two or more identical copies from a given system. Then, one could increase the amount 
of information by performing measurements on the reproduced copies. Such a strategy fails be- 
cause of the no-cloning theorem 0], which is a simple consequence of the quantum superposition 
principle. 

Nevertheless, no principles of quantum mechanics prevent us from characterizing the quan- 
tum state of an ensemble of physical systems. By repeating the measurement many times on 
individual copies, we may arrive at reliable information on the properties of the ensemble. The 
aim of quantum state measurement can be formulated in technological terms: suppose we have 
a machine producing identically prepared copies of a quantum system. When delivering this 
output for some application, we should be able to provide its specification, which on the most 
complete level means a full characterization of the quantum state. Such a problem, apart from 
interesting fundamental aspects, is currently of practical interest in many areas of science. This 
is due to extensive studies devoted presently to preparation, manipulation, and control of quan- 
tum systems. The motivation of this research is to overcome current technological limitations by 
exploiting fully possibilities offered by quantum mechanics. Let us mention just few examples. 
The yield of chemical reactions can be increased by controlling the quantum state of reactants 
[^. Application of so-called squeezed states of light improves precision of interferometric mea- 
surements, which can be used to enhance sensitivity of gravitational wave detectors [^. Coherent 
preparation and manipulation of entangled multiparticle systems can be used to solve computa- 
tional problems intractable by classical computers [^. An important matter in developing these 
and other technologies is the possibility to gain extensive and reliable information on the state 
of quantum systems. The ultimate tool for this purpose is the measurement of the complete 
quantum state. 

When characterizing ensembles, we usually need to take into account the possibility of sta- 
tistical fluctuations, and to describe their state using the more general concept of the density 
matrix rather that a wave function. From the formal point of view, the task of characterizing 
the quantum state can be accomplished by a set of appropriately chosen measurements, which 
would extract unambiguous information on all the elements of the density matrix. The crucial 
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question, however, is: how to do this in practice? The challenge of quantum state measurement 
has several important aspects. The first one is design and realisation of measurement schemes 
that yield complete characterization of the quantum state. Further, there is a nontrivial task of 
extracting precise information on the quantum state from data collected in a realistic, imperfect 
experimental setup. Finally, we have said above that expectation values of a sufficiently large 
set of observables contain complete information on the density matrix. The problem is that true 
expectation values are obtained only in the limit of the infinite number of measurements. In a 
real laboratory we always deal vjith finite ensembles, and this is our source of information on the 
quantum state, which we would like to use as efficiently as possible. 

In quantum optics and related fields, several examples of simple quantum systems have been 
thoroughly studied, including a single light mode, a trapped ion, and a diatomic molecule. A lot 
of interest has been paid to detection of subtle quantum statistical effects. It was therefore natural 
that the domain of quantum state measurement has grown mainly on the ground of advances in 
quantum optics. In 1993, the group of Michael Raymer at the University of Oregon demonstrated 
complete experimental characterization of the quantum state of a single light mode by means of 
optical homodyne tomography [^ . This seminal experiment was followed by extensive research 
in the field of quantum state measurement. Over past several years, we have witnessed a series 
of beautiful experiments with various quantum systems [§] . The vibrational state of a diatomic 
molecule has been reconstructed from measurements of the time-dependent fluorescence spec- 
trum |Q|. Optical homodyne tomography has been applied to reconstruct a whole gallery of 
squeezed states of light [|TOp. The motional state of a trapped ion has been characterized using a 
very sophisticated technique based on the monitoring of the fluorescence [|Tl|]. The tomographic 
method has been used to measure the transverse motional state of an atomic beam [|T2|]. All 
these experiments were tightly connected with the rapid theoretical development of the domain 
of quantum state measurement. Numerous measurement schemes have been proposed and anal- 
ysed in detail. In particular, the role of statistical uncertainty has been discussed, and various 
approaches to reconstructing quantum state representations from experimental data have been 
described. 

The subject of this thesis is the measurement of the quantum state in the phase space. The 
concept of the phase space provides a bridge between the quantum mechanical formalism and 
classical physics. Predictions of quantum mechanics have essentially statistical character. As a 
rule, one can predict only probabilities of obtaining specific outcomes of a measurement. Such 
a situation can be encountered in classical mechanics as well. For example, if we deal with 
an ensemble of classical systems, properties of a single copy can be defined only in statistical 
terms. The state of the ensemble is characterized by a phase space distribution, which describes 
the probability of occupying a given volume element by the system. One may wonder whether 
this intuitive picture of fluctuations has its counterpart in quantum physics. The answer to this 
question is not straightforward. It is possible to convert the quantum mechanical formalism 
into a form which resembles a classical statistical theory. The first phase space representation 
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of the quantum state was introduced in 1932 by Wigner [jT3p. However, such a phase space 



representation is not unique: noncommutativity of quantum observables leads to abundance of 
quantum analogs of the phase space distribution, and none of them captures all the properties of 
the classical object Q]. This is a manifestation of the fact that quantum mechanics is essentially 
different from a classical theory. Nevertheless, quantum phase space quasidistributions contain 
complete information on the quantum state. The family of quasidistribution functions provides 
a convenient framework for studying many quantum optical problems. It is also a useful tool in 
visualising quantum coherence and interference phenomena. 

For a long time, quantum quasidistribution functions have been considered mainly as a quite 
odd theoretical concept rather than a quantity which can be measured in a feasible experimental 
scheme. This perspective changed completely with the demonstration of optical homodyne to- 
mography, which brought quasidistributions, in particular the Wigner function, to the realm of 
a physical laboratory. Optical homodyne tomography is based on the observation that marginal 
distributions of the Wigner function of a light mode can be measured by means of homodyne 
detection. The inverse problem, i.e. the retrieval of the Wigner function from its projections, is 
similar to the procedure used in medical imaging, where the spatial distribution of the tissue is 
reconstructed from absorption measured across the body. Practical implementation of the recon- 
struction algorithm is a rather complex and delicate matter: the back-projection transformation 
is singular, and its application to experimental data has to be accompanied by a special filtering 
procedure. Demonstration of optical homodyne tomography was a successful combination of a 
precise quantum optical measurement with sophisticated data processing. 

In this work we develop a novel, entirely different approach to measuring quasidistribution 
functions of light. We exploit the fact that the value of a quasidistribution at a given point of the 
phase space is itself a well defined quantum observable. Motivated by this representation, we 
propose and demonstrate an optical scheme for measuring directly quasidistribution functions. 
This method, based on photon counting, avoids the detour via complex numerical reconstruction 
algorithm. The basic elements of our measurement scheme are very simple. The light mode 
whose quantum state we want to measure is interfered with an auxiliary coherent probe field, and 
a photon counting detector is used to measure the photocount statistics of the superposed fields. 
We show that a simple arithmetic operation performed on the measured photocount statistics 
yields directly the value of the quasidistribution at a point defined by the amplitude and the phase 
of the probe field. By changing these two parameters of the probe field, we may scan the complete 
phase space, and obtain the full representation of the quantum state of the measured light mode. 
We demonstrate an experimental realisation of this scheme, and present measurements of the 
Wigner function for several quantum states of light. The experimental part of this thesis has 
been performed in Division of Optics, Institute of Experimental Physics, Warsaw University, in 
collaboration with Prof. Czestaw Radzewicz. 

We shall study here in detail various aspects of the direct scheme for measuring quasidistri- 
bution functions. On the practical side, there is a question about the role of typical experimental 
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imperfections. We shall analyse how the result of the measurement is affected by such factors 
as non-unit detection efficiency and imperfect interference visibility. We shall also provide esti- 
mates for statistical error, which are necessary to design an accurate experiment, and to specify 
confidence of the experimental outcome. These theoretical results will be an important tool 
for quantitative analysis of the performed experiment. In addition, our discussion of practical 
aspects has interesting consequences in the recently disputed problem of compensating for de- 
tector losses in photodetection measurements. Over past several years, there were conflicting 
claims concerning the possibility of removing deleterious effects of imperfect detection by ap- 
propriate numerical processing of experimental data [[1^, |T^ [1^]. Our measurement scheme 



provides a testing ground for this problem. We will show that in general no compensation for 
detection losses is possible, unless some a priori knowledge about the measured quantum state 
is given. Discussion of this problem reveals the fundamental role of statistical uncertainty in 
realistic quantum measurements, which results from the fact that in a laboratory we always deal 
with finite ensembles. 

An attractive feature of the presented approach to measuring quasidistribution functions of 
a single light mode is the direct link between the measured observable and the quantum state 
representation. One may wonder whether this approach can be applied in other situations. We 
shall describe here a generalization of the measurement scheme to multimode radiation. We shall 
demonstrate that multimode quasidistribution functions are also directly related to the photocount 
statistics, and that they can be determined in an equally simple way. Although our interest in this 
thesis will be confined to detection of optical radiation, it should be noted that the idea underlying 
our measurement scheme has proven to be fruitful in the measurement of the vibrational state 
of a trapped ion [pT|]. It has also motivated measurement schemes for a cavity mode [ p^ and a 
diatomic molecule [|T9p . 

This thesis is organized as follows. In Chap. ^ we review phase space representations of the 
quantum state. Starting from the definition of the Wigner function, we show how this distribution 
can be unified with other phase space representations, and we discuss properties of generalized 
s-ordered quasidistribution functions. In Chap. ^ we review briefly previous work on measuring 
the quantum state of light. We present two techniques which have been realized in experiments: 
optical homodyne tomography and balanced homodyne detection. Next, in Chap. ^ we intro- 
duce the direct method for measuring quasidistribution functions of light. We present the phase 
space picture of the measurement, and we develop the multimode theory of the scheme. Various 
practical aspects of the proposed measurement are discussed in Chap. ^, including the effects 
of imperfect detection, and the possibility of compensation for detector losses. In Chap. ^ we 
present experimental realization of the proposed scheme, and demonstrate the direct measure- 
ment of the Wigner function of a single light mode. The issue of statistical uncertainty in pho- 
todetection measurements is discussed from a more general point of view in Chap. [7[ We show 
that the statistical noise sometimes limits available information on the quantum state. Finally, 
Chap. |] concludes the thesis. Major part of original results presented in this thesis has been 
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Chapter 2 

Phase space representations of quantum 
state 



In the standard formulation of quantum mechanics, the quantum state is characterized by a vector 
from the Hilbert space describing the physical system. The state vector is related to measurable 
quantities by evaluating expectation values with operators which represent observables. This 
formalism is very far from a classical, intuitive picture of statistical fluctuations. Nevertheless, 
there is a possibility to transform the standard quantum mechanical formalism into the form 
which resembles a classical statistical theory. Such a representation is particularly useful in 
investigating the classical limit of quantum mechanics. The fundamental role in this approach 
is played by quasidistribution functions, which can be regarded as quantum analogs of a phase 
space probability distribution. However, due to noncommutativity of quantum observables, the 
phase space representation of the quantum state is not unique, and it is not possible to have in 
quantum mechanics a phase space distribution that has all the properties of the classical one. 

2.1 Wigner function 

In 1932, Eugene Wigner [ [T3| ] introduced a quantum analog of the classical phase space proba- 
bility distribution. For a particle travelling along one dimension, the Wigner function is related 
to the wave function ip{x) through the formula: 

W{q,p) = -^ /dxV^*(g + x/2)e'^"/'^7/'(g-x/2), (2.1) 

Znn J 

and it completely characterizes the quantum state. The integral of W{q,p) over q and p is one, 
which follows from the normalization of the wave function. Expectation values of quantum 
observables can be obtained from the Wigner function by integrating it with appropriate Wigner- 
Weyl expressions representing these observables [[Hp. Furthermore, marginals of the Wigner 
function yield quantum mechanical distributions for the position and the momentum. However, 
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the Wigner function has one property which manifests that quantum mechanics is distinct from a 
classical statistical theory: the Wigner function can take negative values. We shall see later that 
this property is closely related to quantum interference phenomena. 

Difficulties with defining the quantum phase space distribution have their origin in the non- 
commutativity of quantum observables. As the position and momentum operators do not com- 
mute, we cannot introduce a joint distribution of these two observables. This problem is closely 
related to the issue of the ordering of observables, which appears when passing from classical to 
quantum mechanics. For example, the classical expression qp has the following nonequivalent 
quantum counterparts: qp, pq, or |(gp + pq). The Wigner function corresponds to a specific, 
symmetric ordering of the position and momentum operators, called the Weyl ordering [PD[]. We 
will now transform Eq. ( [2.ip to the form which shows explicitly relation between the Wigner 
function and the symmetric ordering of the position and momentum operators. For this purpose, 
let us introduce an additional delta function and represent it in an integral form: 

W{q,p) = ^JdxJdye^P-/%y-q-x/2)r{ymy-^) 

= j^ jdx j dy jdk e^^^/V(^-^-^'/2)'=/V*(y)V^(t/ - x) 

= — ^ /"da; /"dfce'^P"-^'?)/'' /"dye''="/2;i^*(y)e''=(^-")/'*V'(l/-a;). (2.2) 
In the last expression, the integral over y can be written as the quantum expectation value: 

^^gifcx/2?i^*(^^)gifc(j/-x)/a^(^ - x) = (V^|e*('=^~^^)/^|V^). (2.3) 

This quantity is a function of two real parameters k and x. By differentiating over k and x 
we may obtain moments of the position and momentum operators. These moments are ordered 
symmetrically in q and p, which follows from the form of the exponent in Eq. (|2.3|). The function 
^^|g«(fcg-px)/?i|^^ is called the Wigner- Weyl ordered characteristic function for the position and 
the momentum. Coming back to Eq. (|2.2D, we finally arrive at the formula 



W{q,p) = — ^ f dx I dke'^P''-^'i^/^{i)\e'^^^-^''^'^\i)) (2.4) 

which shows that the Wigner function is the Fourier transform of the symmetrically ordered 
characteristic function for the position and the momentum. Eq. ( |2.4[ ) can be used to evaluate the 
Wigner function corresponding to a mixed state described by the density matrix q. In such a 
case, we have to replace i^'t\)\e^'^^^~^^^l^\'\\)) by Tr(^e'(*^'^~^^'^/^). Equivalently, the Wigner function 
of a mixed state can be obtained from a weighted sum of the Wigner functions describing the 
pure components of the mixed state. 

For a harmonic oscillator, it is convenient to introduce a pair of dimensionless annihilation 
and creation operators, defined by the equations: 

a=^{\-^q^%\lr^p), a^ = ^\-^q-iXh-^p) (2.5) 

V 2 v2 
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where A is a natural length scale defined by the mass and the frequency of the oscillator. Anal- 
ogously, the two real parameters of the Wigner function can be combined into a single complex 
argument a = {X^^q + i\fV^p)/\/2. The Wigner function in this parameterization is given by 



1 



Wia) = 4 / d\e'^*°-^"*(e«'^^-^*'^) (2.6) 

~r-2 I 



where the integration is performed over the whole complex plane and the angular brackets (...) 
denote the quantum expectation value. Let us note, that the normalization constant in Eq. ( |2.6| ) 
has changed compared to Eq. (|2.4|). This is because the integration measure over the phase space 
is now equal to d^a = dq dp/2h. 

The physical system which we shall describe in the phase space representation, is optical ra- 
diation. In the standard procedure of quantization, the electromagnetic field is decomposed into 
a set of independent modes. Each of these modes is characterized by a pair of creation and anni- 
hilation operators, which satisfy bosonic commutation relations for a harmonic oscillator. When 
only one of the modes is excited, we may describe its quantum state using the Wigner function 
defined in Eq. (|2.6|). In the classical limit, the parameter a characterizes the complex amplitude 
of the field, expressed in dimensionless units. We may also use the single-mode description if 
our measuring apparatus is sensitive only to a selected mode of the detected radiation. 

2.2 Quasidistribution functions 



As it is clearly seen from Eq. (p^, the Wigner function corresponds to the characteristic func- 
tion with the symmetric ordering of the creation and annihilation operators. In principle, we 
could consider also other orderings, for example normal or antinormal. In the normal ordering 
all creation operators are placed before annihilation operators, and vice versa for the antinor- 



mal ordering. Thus, we could think of replacing the quantum expectation value in Eq. ( [2.6D 
by the normally ordered characteristic function (e''° e^''*"), or by the antinormally ordered one 
(e~^*"e^" ). These and other possibilities can be written jointly in a very elegant way by introduc- 
ing an exponential factor, which defines the ordering of the creation and annihilation operators. 
This idea leads to the concept of more general s-parameterized quasiprobability distributions. 
The one-parameter family of quasidistribution functions is given by the following formula [P|]: 



,s|c|2/2+C*"-C"* /^Cat-ca 



W{a;s) = ^ f d^C e^l^l^/^+^'°-^°' (^e^'^' '^''^^ . (2.7) 

The real parameter s is associated with the ordering of the field bosonic operators through the 
exponential factor exp(s|Cp/2). In particular, it can which can be easily checked using the 
Baker-Campbell-Haussdorf formula that three values s = 1,0, and —1 generate the normal, 
symmetric and antinormal ordering, respectively. 
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The definition given by Eq. ( |2.7[ ) unifies the Wigner function with other, independently devel- 
oped quantum analogs of a phase space distribution. For instance, normal ordering corresponds 



to the so-called P function, introduced by Glauber [ p2| ] and Sudarshan [|23|]. This function serves 



as a weight function in the diagonal coherent state representation for the density matrix g: 

d^aP{a)\a){a\. (2.8) 



On the other hand, antinormal ordering yields the distribution known as the Husimi [ |24| ] or Q 



function [ P5| , |2q ], which is given by the diagonal elements of the density matrix in the coherent 



state basis: 

Q{a) = -{a\g\a). (2.9) 

Properties of various s-parameterized quasidistribution functions are quite different. This 
can be seen using the three examples of the P function, the Wigner function, and the Q function. 
The P function is highly singular for nonclassical states of light. For example, it is given by 
derivatives of the delta function for eigenstates of the photon number operator a^a. The Wigner 
function is well behaved for all states, but it may take negative values. Finally the Q function 
is always positive definite, which follows directly from Eq. ([2.9[). The fact that quasidistribution 
functions with lower ordering are more regular reflects a general relation linking any two differ- 
ently ordered quasidistributions via convolution with a Gaussian function in the complex phase 
space: 

where s > s'. Thus the lower the ordering, the smoother the quasidistribution is, and fine details 
of the function are not easily visible. The Gaussian exponent appearing in the above equation 
can be formally regarded as a propagator for the diffusion equation, with the ordering parameter 
playing the role of the time. Following this analogy, we may write a differential equation for 
quasidistribution functions corresponding to a given quantum state: 

Let us note that the above equation differs from the standard diffusion equation by the minus 
sign. This difference originates from the fact, that the "diffusion" of quasidistributions follows 
in the direction of decreasing s. 

In our calculations a normally ordered representation of the quasidistribution functions will 
be very useful. Introducing normal ordering of the creation and annihilation operators in Eq. ( |2.7[ ) 
allows to perform the integral explicitly, which yields: 



W{a; s) = " ^^ ( : exp ( -^^^(a^ - a*){a - a) ] : ) . (2.12) 
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Thus the s-ordered quasidistribution function at a complex phase space point a is given by the 
expectation value of the operator 

W(a; s) = — : exp I (a^ - a*)(d - a)] : (2.13) 

7r(l — s) \ 1 — s J 

Using the operator identity p7| ] 

: exp[(e*'' — l)-0'{;] : = exp{i(^v^v) (2.14) 



valid for an arbitrary bo sonic annihilation operator v, we may transform Eq. ( |2.13D to the follow- 
ing expression: 

2 /s + 1^ (at-a*)(a-a) 



Wia;s) = — . (2.15) 

7r(l — s) \s — 1/ 

The operator appearing in the exponent is the displaced photon number operator n = a) a. Using 
the standard displacement operator D{a) = exp(a;d^ — a*d), we may write W{a] s) as [ p8| , ^ : 

= -7, vE^ D{a)\n){n\D\a). (2.16) 

7r(l - s) ^\s-lj 

^ ' n=0 ^ ' 

The last form is simply the spectral decomposition of Wifit; s). The eigenvectors are displaced 
Fock states D{oi)\n), and the corresponding eigenvalues are [(s + l)/(s — 1)]" times the front 
normalization factor 2/7r(l — s). It is instructive to see, how the properties of the quasidistri- 
butions are reflected by the spectrum of iy(a; s). First, let us note that for s ^ 1, the factor 
(s + l)/(s — 1) is divergent; this corresponds to the singular character of the P representation. 
For < s < 1 the set of eigenvalues is unbounded; therefore, the corresponding quasidistribu- 
tions also may exhibit singular behaviour. The operator Wia; s) becomes bounded for s < 0, 
and the highest value, i.e. s = 0, corresponds to the Wigner function. Even when W{a\ s) is 
bounded, its eigenvalues corresponding to odd ns can be negative. The highest value of s for 
which all the eigenvalues are nonnegative is s = — 1, which corresponds to the Q function. 

2.3 Quantum interference in phase space 

We will now discuss, using a simple example, how quantum interference phenomena are visu- 
alised in the phase space representation. A quantum analog of a classical field with well defined 
amplitude and phase is the coherent state |«o)5 defined as an eigenstate of the annihilation op- 
erator a\oLQi = ao\ao). This equivalence originates from the fact, that full quantum theory of 
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photodetection gives for coherent states the same predictions as semiclassical theory with quan- 
tized detector and classical electromagnetic fields [^. Coherent states are represented in the 
phase space by Gaussians 



p|/l"o>( 



TT 1 



exp 



1 -s' 



\a 



ttol 



(2.17) 



Quantum mechanics allows one to combine two such classical-like state, for example \ao) 
and I — «o) into a coherent superposition 



v/2(l + e-2|"oP) 



(|ao) + I -tto))- 



(2.18) 



States of this type illustrate quantum coherence and interference between classical-like com- 
ponents, and are often called quantum optical Schrodinger cats [|T]]. In contrast to coherent 



states, they exhibit a variety of nonclassical properties [p2|]. The quasidistribution function of the 
superposition I?/)) is given by the formula 



iyl^>(a; 



1 



7r(l-s)(l + e-2KP) 
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exp 



1-s 



|q; — aol 



+ exp 



+2 exp 



1-s 
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1-s' 



« + ao 



«o| exp 



1-s' 



2 \ f Alm(aoa*) 
\a\ cos 



V 1-s 



. (2.19) 



The first two terms in the square brackets describe the two coherent components. The last term 
results from quantum interference between these components. It contains an oscillating factor 
cos[4Im(aoa;*)/(l — s)]. It is seen that the frequency of the oscillations grows with the dis- 
tance between the coherent components. The envelope of this oscillating term is defined by the 
Gaussian exp[— 2|ap/(l — s)], which is centered exactly half way between the interfering states. 



Fig. ^^ shows quasidistributions plotted for three different values of the ordering parameter 
s. The Wigner function contains an oscillating component originating from the interference be- 
tween the coherent states. This component is much smaller for s = —0.1 and it is completely 
smeared out in the Q function, which can hardly be distinguished from that of a statistical mix- 
ture of two coherent states. This is because the whole interference term is multiplied by the 
factor exp[2s|Q;oP/(l — s)], which quickly tends to zero with decreasing s. Let us note that the 
larger is the distance between the components, the faster this factor vanishes. The decay of the 
interference component can be formally viewed as a result of diffusion, described by Eq. ( [2. lip . 
Decreasing the ordering parameter s makes the whole quasidistribution blurred, and this effect 
is particularly deleterious to the quickly oscillating pattern. 
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Figure 2.1: Quasidistributions representing the Schrodinger cat state for oq = 3i, depicted for the ordering 
parameters s = 0, —0.1, and —1. 

2.4 Deconvolution 



We have seen, using the example of the Schrodinger cat state, that signatures of quantum inter- 
ference can be visible better in quasidistribution functions with higher ordering. Thus, what is 
interesting, is the inversion of Eq. ( [2.10D , i.e. deconvolution of a lower-ordered quasidistribution 
function. This task is quite difficult. Let us first note that in general the integral in Eq. ( |2.10D fails 
to converge if we take s < s'. Instead, we may use the Fourier transforms of the quasidistribution 
functions 



W{C,s) 



d^pe 



C0*-C 



^W{(3;^ 



s\C\V2,ra^-C& 



)• 



(2.20) 



Transition to a higher ordered quasidistribution consists now simply in multiplication by an ex- 
ponent: 



W{(- s') 



.{s'-sm- 



^'W{C; 



(2.21) 



and evaluation of the inverse Fourier transform. The complete expression of 1^(0;; s') in terms 
of a lower ordered quasidistribution has the form: 



W{a;s') 



vr^ 



^2^^{s'-s)\(\y2+Ca-Ca* 



d^f3e^^*-^'^W{(3;s). 



(2.22) 



Anticipating for a moment the connection of the quasidistributions with experiment, let us sup- 
pose that we are given an experimentally determined quasidistribution W{P] s), and that we are 
trying to apply the deconvolution procedure described by Eq. ( |2.22| ). Usually, values of W{P; s) 
will be affected by errors originating from statistical uncertainty and various experimental imper- 
fections. These errors make the deconvolution a very delicate matter. The crucial problem is that 
the Fourier transform W{C; s) has to be multiplied by an exploding factor e*^"* ^''^'''l /^. Experi- 
mental errors of W{P; s) can generate long, slowly decaying high-frequency components in its 
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Fourier transform. Multiplication by an exploding exponent enormously amplifies contribution 
of these fluctuations, which leads to huge errors of the reconstructed W{a; s). Therefore, decon- 
volution of experimentally determined quasidistributions according to Eq. ( |2.22D is practically 
impossible. 

2.5 Multimode quasidistributions 

The concept of quasidistribution functions can be generalized in a straightforward manner to 
multimode radiation. In analogy to Eq. (p77|), we need to take the symmetrically ordered multi- 
mode characteristic function, and to evaluate its Fourier transform with an appropriately chosen 
Gaussian factor which defines the ordering: 

W{ai,...,aM;s) 

^ / dCi . . . dCM exp ( XI 1 1^*1^ + ^i^i ~ ^i^i ) ( ^^P ( 5Z ^^^i ~ ^i' 



^2M 



(2.23) 

Introducing normal ordering allows one to perform the integrals, which yields an explicit nor- 
mally ordered representation: 

W{ai,...,aM;s) = (———] {■.expi-Y^^{al-(^i)io^i-ai)\:\. (2.24) 
Using Eq. ([2.14p, we may represent the quasidistribution functions as: 



^■^ / /S + 1\ ^fiM-»m^-'^^)' 



The expression X]i=i(^l ~ a*)(oj — ««) appearing in the exponent is simply the phase space 
displaced operator of the total number of photons. In analogy to the single-mode case, we may 
write the multimode quasidistributions as an expectation value of the operator involving the 
multimode displacement operator 

D{{ai}) = exp ^aiO,] - a*ai j , (2.26) 



.4 = 1 



and the total photon number operator, defined as 



M 



N = J2^i (2.27) 



4 = 1 
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where hi = dldi. The explicit expressions are: 



2 \ ^ V ^ . . . . / 2N 



W{ai,...,aM;s) = [—- { D{{ai}) : exp 

7r(l — s) J \ \ 1 — s 

Ml / I 1 \ * 



vrfl — s) J \ \s — 1 




2.6 Quasidistribution functionals 

The representation given in Eq. (|2.2g|) suggests generalization of the multimode quasidistribution 



functions to the form independent of the specific decomposition into modes. Such generalized 
quasidistributions are functionals of the electromagnetic field. Instead of using a finite set of 
annihilation and creation operators, we will now deal with the full description of the electromag- 
netic field, involving the operator fields E(r, t) and H(r, t). In order to simplify the notation, we 
shall fix the time t, and omit it in the subsequent formulae. The definition of quasidistribution 
functionals involves two operators: the coherent displacement operator V, and the operator of 
the total number of photons J\f. The action of the displacement operator is straightforward: it 
adds a classical amplitude to the field operators according to the formula 

P[E(r),H(r)]E(r)I'T[E(r),H(r)] = E(r) - E(r) 
:D[E(r),H(r)]H(r)I'"f[E(r),H(r)] = H(r)-H(r). (2.29) 

In order to find an explicit formula for quasidistribution functionals, we need to express the total 
photon number operator J\f in terms of the electric and magnetic field. We shall start from the 
standard decomposition of the electromagnetic field into plane waves with periodic boundary 
conditions in a box of the volume V: 



EW = ^Ev|^^''^(«''^^^''^-«L^"^''^) (2.30) 

la 



H(r) = -— EV^^'^x^(^'^^^'"-^l^"^'")- (2.31) 

c/io ^ V 2eoF k/ 

t(j 

Here the indices / and a label respectively the wave vectors k^ and the polarizations e^^r, and 
uji = c\ki\ is the frequency of an /th mode. Our goal is to represent the sum 



■A^=$^«l>- (2-32) 



ho- 



using E(r) and H(r). For this purpose we shall take Fourier transforms of these fields: 



d^rEWe-^^'"- = ^^|^J2i^^^^l^-^-l^^~J (2-33) 
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/d^rH(r)e-^^- = -m/ ^^ E f^'^ >< ^«'^ " ^-'^ x ^«-^^V (2.34) 

Here on the right-hand sides we have used the fact that k ^ = — k;. The product of the Fourier 
transforms taken for k; and k_i can be expressed as: 

(2.35) 
and 

d=^r /d-\'H(r)H(r')e-*'^'('-'') 



fiioiV -^-^ f k; „ k_/ ^i \ f ki „| k_i „ \ 

= -7: y2{(^aa'aiaa]^, + ei^d^i^icii^a^i^: + e-ia^ia'a^^i^a]^, + <5aa'al;^a_;^0- (2.36) 

crcr' 

We shall now add the expressions for the electric and magnetic fields multiplied by the factors 
eQ/2huJiV and fio/2JTiuiV respectively. This yields: 

fdh fdh' (|E(r)E(rO + ^H(r)H(rO) "'^V ^ = ^ ^^"'-"''- + ^-i.^-i^)- (2-37) 
J J I ^ 

This formula is close to the standard expression for the energy of the electromagnetic field. 
Indeed, we could obtain it via multiplication of both the sides by huj, and summation over /. 
However, we are now interested in a different quantity, namely the total number of photons, and 
we need to perform the summation with the factor hw in the denominator of the left hand side. 
In this way we obtain: 

E l(«''^«l + «f>.) = /d'r/d'r' (|E(r)E(rO + |^H(r)H(r')) K{v - v'). (2.38) 
la- 
in the second term of the left-hand side we have changed the summation index —l^l. The 
integral kernel K(r — r') appearing on the right-hand side is given by: 

We shall evaluate it in the continuous limit, when the sum over / can be replaced by a three- 
dimensional integral over the wave vector k. In this limit, there occurs a singularity at r = 0, 
which is a result of the slowly decaying integrand with large k. We shall regularize the integral 
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by introducing the upper cut-off k^^x for the wave number. Physically, this means that we do 
not take into account photons with energy larger than hckmax- The regularized kernel can be 
evaluated in a straightforward manner: 

Kir) -- 
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/d^k^"^'^ 


1 rkiimx rT^ 

/ c\hh / rl?9 "in ?9r"*'^l'"''=°^'' 


{2'Kyhc J 
1 1- 
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{27iyhc Jo Jo 


271^ he 


r2 





(2.40) 

It is easily seen that truncation of the wave vector magnitude has removed singularity of the 
kernel K(r) occurring at r = 0. 

On the left-hand side of Eq. ( [2.38D , we have a symmetrically ordered product of the creation 



and annihilation operators ^{aiaaj^ + aj^ai^)- In order to obtain the total photon number oper- 
ator, we need to introduce the normal ordering of the right-hand side of Eq. ( [2.38p . Using this 
expression, we can easily define quasidistribution functionals of the electromagnetic field. There 
is a small difficulty arising from the fact that we now deal with the infinite number of degrees 
of freedom. In Eq. ( |2.24j ), we cannot pass to infinity with the number of modes in the normal- 



ization prefactor [2/7r(l — s)]^^ . We shall solve this difficulty by absorbing the normalization 
prefactor into the functional integration measure over the fields E and H. Thus, we define the 
quasidistribution functional as: 

W[E(r),H(r);s] = /p[E(r), H(r)] : exp [-^] : ^^[E(r), H(r)] \ . (2.41) 

This definition can be written explicitly using the electromagnetic field operators with the help 
of the derived expression for the total photon number operator as: 

W[E(r),H(r);s] 

2 f ,3_. f j3_./ r^/_. .J\ /^O 



exp 



1-s 



fd^r fd^r'K{r - r') (-[E(r) - E(r)][E(r') - E(r') 



+ ^[H(r)-H(r)][H(r')-H(r')]) 



(2.42) 



The integral kernel appearing in the above formula is given by Eq. ( |2.39D . 
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Chapter 3 

Homodyne techniques for quantum state 
measurement 



Over the last decade, the domain of quantum state measurement has passed a long way from 
first theoretical proposals to well understood experimental realizations. Complete presentation 
of the current state of this field would require a separate book, encompassing a wide range of 
experimental techniques and concepts of data analysis. In this chapter we shall set the scene for 
further parts of the thesis by describing briefly earlier works on measuring the quantum state of 
light. We shall restrict our attention to detection of optical radiation, and describe two techniques 
which have been successfully realized in experiments: double homodyne detection and optical 
homodyne tomography. 

Double homodyne detection allows one to measure the Q function of a light mode. It was 



demonstrated in 1986 by Walker and CaroU p3|]. Their experiment had as a main purpose the 



demonstration of a homodyne measurement near the quantum noise limit, and it later attracted 
attention as a complete characterization of the quantum state. As we discussed in the previous 
chapter, the Q function is a positive definite distribution, and it exhibits only faint traces of 
quantum interference. Optical homodyne tomography was realized first by Smithey et al. in 
1993 [0]. This technique is capable of measuring the Wigner function. Apparently, this fact 
added extra excitement to the development of homodyne tomography, as the Wigner function 
is a nonclassical distribution function which may take negative values resulting from quantum 
interference. 

Both these techniques are based on the same experimental apparatus, namely the balanced 
homodyne detector. This device provides information on phase-sensitive properties of light. In 
the quantum mechanical formalism, it performs the measurement of a family of observables 
called quadratures. We shall start this chapter with a description of the balanced homodyne 
detector in Sec. p7T[ The double homodyne detection scheme is discussed in Sec. ^7^ . We show 
that this scheme can be used to measure two noncommuting observables at the cost of introducing 



extra noise to the measurement. Sec. ^3] is devoted to optical homodyne tomography. It describes 
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the physical principle of the method, as well as mathematical transformations involved in the 
processing of experimental data. 

3.1 Balanced homodyne detector 

Standard photodetection is insensitive to phase properties of optical radiation. This is because the 
observed signal depends only on the operator of the number of photons n = a) a. Nevertheless, 
we may use a photodetector to measure phase-dependent quantities by superposing the measured 
beam with an auxiliary coherent field using a beam splitter. The auxiliary field has the name of 
the local oscillator. When measuring such a superposition, the signal from the photodetector will 
be described by an expression involving terms linear in a and a) . Thus, it carries information on 
the phase properties of the measured field. This is the basic idea of homodyne detection. 

In quantum optics, homodyne detection has played an important role in investigating the 
squeezed states of light. These states exhibit interesting noise properties in certain phase-de- 
pendent observables. More specifically, for a single light mode we may introduce a family of 
quadrature observables dependent on the phase 9: 



e'^a) + e 



(3.1) 



It is easy to check that the commutator of two quadratures corresponding to phases which differ 
by 7r/2 is [xq,xqj^t,/2\ = i. Consequently, variances of these two observables satisfy the un- 
certainty relation in the form Ax0Axg+T,/2 > 1/2. For coherent states, this variance is evenly 
distributed over all the quadratures, and it equals to Axq = 1/v^. Squeezed states are such 
states of the electromagnetic field, which for a certain phase 9 have the variance smaller than the 
coherent state level. These states cannot be described within classical theory of radiation, and 
the squeezing is clearly a non-classical property. Squeezed states can find application in very 
precise interferometric measurements [^ . 

Of course, it is obvious from the definition of quadratures that their measurement requires a 
phase-sensitive technique, such as homodyne detection. In practice, we need to take into account 
various experimental imperfections. One of them is the excess noise of the local oscillator field. 
This noise adds to the observed level of fluctuations, and it may mask subtle quantum effects 
related to squeezing. Fortunately, there is a possibility to subtract the local oscillator noise by 
using the so-called balanced scheme. In the balanced homodyne detection scheme, depicted 
in Fig. PTT] , the signal field is superposed using a 50:50 beam splitter with the local oscillator. 
Two photodetectors monitor the output ports of the beam splitter, and the recorded signal is the 
difference of the detector photocurrents. In this way, we cancel the effect of local oscillator noise 
when measuring the variance of the difference photocurrent. 

Moreover, theoretical description of the balanced homodyne detector shows that this setup, 
in an idealized limit, is the optical realization of the quantum measurement of the quadrature op- 
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50:50 



Figure 3.1: The balanced homodyne setup. The signal field, described by a density matrix p, is combined 
with a coherent local oscillator |/3)lo- The two outgoing fields are measured using photodetectors. The 
difference of their counts is the statistical data recorded in the experiment. 



erator xq. The idealization is based on two assumptions: the unit efficiency of the photodetectors 
and the classical limit of the local oscillator. The latter condition can be easily satisfied in an 
experiment. Also, efficiency of photodetectors used in a homodyne setup can be close to 100%. 
If we rely on these two assumptions, theoretical analysis of the balanced scheme becomes quite 



compact. More detailed studies can be found in Ref. [ |34| , |35|, |3q]. 

Let us describe the signal field with the annihilation operator d, and the local oscillator field 
with b. These two fields, superposed on a 50:50 beam splitter, yield two outgoing modes. In 



general, combination of the modes at a beam splitter is given by an SU(2) transformation [ pTQ 
For a 50:50 beam splitter, we may simplify this transformation to the matrix 



1 

72 



(3.2) 



where c and d are the annihilation operators of the outgoing fields. We assume that the local 
oscillator is in a coherent state \P)lo^ and the quantum state of the mode d is given by the 
density matrix g. 

The quantity we are interested in is the difference of photocurrents generated by the detectors 
monitoring the modes c and d. On the microscopic level, these photocurrents consist of a discrete 
number of electrons rii and ^2. In a real experiment, this discreteness is not observed due to the 
large average number of the generated photoelectrons. The observable measured in balanced 
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homodyne detection is the difference of the electron number AA^ = rii — n2. The probability 
p(AN) of obtaining a specific value for AA^ can be easily derived using the standard theory of 
photoelectric detection. It is given by the expression: 

p{AN)= V Tr{g®\f3){(3\Lo--e-'''^-^e-'''^-^-^:}. (3.3) 

AM '^l' ^2! 

ni—n2=i\N 

In further calculations, it is more convenient to deal with the generating function for the probabil- 
ity distribution p(AA^). The generating function is obtained by evaluating the Fourier transform: 

00 

AN=—oo 

= Tr{g ® \P){P\lo ■ exp[(e^« - l)c^c + (e-^« - l)d^ : }. (3.4) 

In this way, we managed to get rid of the troublesome constrained sum with the condition rii — 
ri2 = AA^. We can now remove the normal ordering symbol by making use of the operator 
identity given in Eq. ( |2.14D . This yields: 

Z(0 = Tr{g^\(3){(3\Loe'^'-'''-'^''^^} = Tr{Q ^ \f3){f3\Loe'^^'''+'''^}. (3.5) 

When the local oscillator is in a strong coherent state, the bosonic operators b, W in the exponent 
gr5(a b+ah ) ^^^ ^^ replaced by c-numbers /5, /?*. In this regime, it is also convenient to rescale 
the difference photocurrent AA^, which grows as the first power of the local oscillator amplitude. 
Dividing AA^ by \(3\, we obtain a quantity which is independent of the magnitude |/3| in the 
regime of the classical local oscillator. We shall introduce an extra factor of 1/v^, and define 
the homodyne variable as a; = AN / \/2\l3\. This variable can be treated as a continuous one, as 
the local oscillator amplitude is very large. The rescaling of the homodyne variable corresponds 
to changing the parameter of the generating function according to A = ^v^|/3|. In the new 
parameterization, the generating function takes the form: 

Ze{\) = /exp (^{e^'a^ + e-''a)\ \ = (e^^^«), (3.6) 



where 6 is the phase of the local oscillator: (3 = |/5|e*^. We have added here a subscript 6 to the 
generating function Zq{\) to stress that the measured observable depends on the local oscillator 
phase. Using the last form of the Zg{\), we may easily obtain the probability distribution pe{x) 
for the homodyne variable x by evaluating the inverse Fourier transform. Let us note that we 
should now integrate over all real values A because of the introduced rescaling. The inverse 
Fourier transform yields: 

pe{x) = — / dAe-*"-'(e'^"«) = (^(x-xe)) = {\x)ee{x\) . (3.7) 

This expression clearly shows, that balanced homodyne detection is the measurement of the 
quadrature operator xq. The probability of obtaining the result x is given by the projection on the 
corresponding eigenstate of the quadrature operator, defined as X0\x)e = x\x)g. 



CHAPTER 3. HOMODYNE TECHNIQUES. . . 25 

3.2 Double homodyne detection 

In homodyne detection, we have to select the phase of the local oscillator, which defines the mea- 
sured quadrature. Simultaneous measurement of different quadratures is not possible, because 
they correspond to noncommuting observables: it is easy to check that for example [xq, xg+Tr/2] = 
i. However, we may try to circumvent this difficulty by splitting first the input beam on a 50:50 
beam splitter and performing two homodyne measurements on the outgoing fields. This is the 
idea of double homodyne detection. The corresponding setup is shown in Fig. ^J\ The signal 
field, described by an annihilation operator a, is divided using the 50:50 beam splitter BS. The 
two outgoing beams are measured with two separate balanced homodyne detectors. The phases 
of local oscillators can be independently adjusted in each of the arms of the setup, which allows 
one to measure two arbitrary quadratures of the fields leaving the beam splitter BS. For simplic- 
ity, let us choose the two local oscillator phases to be and 7r/2, and to denote the corresponding 
quadratures hy q = xq and p = x^/2. These two quantities commute to the imaginary unit 
[q, p] = i, and they are optical analogs of the position and the momentum operators for a particle. 
In the quantum description of the setup, we need to take into account the vacuum field enter- 
ing through the unused input port of the beam splitter BS dividing the signal field. This vacuum 



field is denoted with the annihilation operator v in Fig. ^^. The quadratures measured at the two 
homodyne detectors are given by the combinations 

Qi = ^{Qa + Qv), P2 = -7k(p^ " P^) (3.8) 

where the indices a and v denote quadrature operators corresponding the signal and the vacuum 
mode respectively. The form of these combinations follows directly from Eq. (|3.2[), describing 
transformation of the field operators at a 50:50 beam splitter. 

The probability distribution p{qi,P2) for the outcomes of the measurement is now defined 
on the two-dimensional space spanned by the variables gi and p2. Analogously to the previous 
section, it will be more convenient to use the generating function Z{Xi, A2), which depends now 
on two parameters Ai and A2: 

^(Ai,A2) = y"dAiy"dA2e*"^'^^+^^^«Mgi,P2). (3.9) 

The generating function describing the joint measurement of the quadratures gi and p2 is given 
by a straightforward generalization of Eq. (|3.6D: 



Z(Ai,A2) = {exp{iXiqi+iX2p2))a,: 



= exp ( --(Af + A^) 1 / exp {-J^i^iQa + ^2pa) j ) • (3.10) 

In the second line we have evaluated explicitly the quantum expectation value over the vacuum 
mode. The joint probability distribution p(qi , P2) can be obtained from the double inverse Fourier 
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Figure 3.2: Double homodyne detection setup. The signal field, denoted by the annihilation operator a, is 
divided using a 50:50 beam spUtter BS. The two outgoing fields fall onto balanced homodyne detectors. 
The local oscillator phases are adjusted such that two conjugate quadratures are measured. In the quantum 
mechanical description of the setup, one has to take into account the vacuum field v entering through the 
unused input port of the beam splitter BS. 



transform of the generating function. We shall rearrange this expression to the form: 

-|ClV2+C*(9i+«P2)-C{9l-iP2)/gC'it-C*a\ 



(27r) 
1 



n 



d\e- 



(3.11) 



where we have substituted C = (zAi — A2)/2. The last expression can be directly related to the 
definition of quasidistribution functions in Eq. (|2.7D, with a = qi + ip2, and s = —I. Thus, 
the joint probability distribution of homodyne events measured in double homodyne detection is 
equal to the Q function of the mode d: 

piQi,P2) = Qa{qi + m)- (3.12) 

One may wonder how this formula changes when we inject an arbitrary state in the second 
input port of the beam splitter BS dividing the signal field. In this case Z(Ai, A2) can be factorized 
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to the product of the symmetrically ordered characteristic functions for the position and the 
momentum: 

Z{Xi, A2) = / exp ( —{XiQa + X2Pa) j) ( exp ( — (Aig„ + X2Pv) j ) ■ (3.13) 

The inverse Fourier transform maps the product of the symmetrically ordered characteristic func- 
tions onto a convolution of the corresponding Wigner functions. After a simple calculation, we 
obtain: 

p{qi,P2) = 2 fdq f dpWa{q,p)W,{V2qi - q, V2p2-p), (3.14) 

where Wa{q^ p) and W^iq^ p) are the Wigner functions describing the quantum state of the fields 
incident on the beam splitter BS. 

The above results illustrates the operational approach to the joint measurement of the position 



and the momentum [ |38| , p9| ]. These two observables do not commute and they cannot be mea- 
sured simultaneously. Nevertheless, we may introduce an auxiliary system, called the "quantum 
ruler", and measure two commuting combinations of positions and momenta. Such a pair of 
combinations has been defined in Eq. (p.8D. These two operational observables can be detected 
simultaneously, and their measurement yields a joint two-dimensional probability distribution of 
two variables which can be related to the position and the momentum. The resulting operational 
phase space distribution is given by a convolution of the Wigner functions of the measured sys- 
tem and the ruler. Double homodyne detection is an optical realisation of this approach, with 
the role of the quantum ruler played by the vacuum field. The vacuum field is described by the 
gaussian Wigner function, and the double homodyne detection yields a smeared Wigner function 
of the signal field, which coincides with the Q function. 



Double homodyne detection has been realized experimentally by Walker and CaroU [|33|]. 
A thorough discussion of this technique can be found in the article by Walker ^^. The same 
experimental scheme has been applied in the operational measurement of the quantum phase 
PT|], and it was shown later that the phase distribution measured in this scheme corresponds to 
the radially integrated Q function 



3.3 Optical homodyne tomography 

In contrast to the Q function, the Wigner function does not have the operational meaning of a 
probability distribution, simply because it may take negative values. Therefore, one cannot de- 
sign an experiment, in which the joint statistics of two real variables would be described by the 
Wigner function. Nevertheless, one-dimensional projections of the Wigner function are posi- 
tive definite. Furthermore, these projections describe quadrature distributions according to the 
formula: 

P0{x) = dyW{xcos6~ysm9,xsm9 + ycos6). (3.15) 
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Figure 3.3: Quadrature distributions pe{x) are one-dimensional projections of the Wigner function 
W{q,p), obtained by integrating it along parallel stripes in the phase space. 



This equation is a generalization of the marginal properties of the Wigner function for the posi- 
tion and the momentum. As we have seen in Sec. P7T| , quadrature distributions can be measured 
by means of a balanced homodyne detector. Thus we may obtain experimentally a family of 
one-dimensional projections of the Wigner function, depicted schematically in Fig. P3| . What 
we would like to do, is to reconstruct from these "shadows" the two-dimensional Wigner func- 
tion. 

A very similar problem is encountered in medical tomography [Q]. By measuring absorption 
of radiation across the body, we can obtain density of the tissue integrated along the direction of 
the measurement. These data are subsequently processed to reconstruct the full density distribu- 
tion using numerical back-projections algorithms. In analogy to this technique, the method of 
reconstructing the Wigner function from homodyne statistics has been called optical homodyne 
tomography. 

The relation between pe{x) and W{q, p) defined in Eq. ( |3.15| ) can be viewed as a transforma- 
tion between two-dimensional functions according to the formula: 



Pe{x) 



dq I dip5{x — qco^O — psm.9)W{q,p), 



(3.16) 



which is known in the field of image processing as the Radon transform. A function W{q.,p) 
of two real variables q and p is transformed into another function pe{x) which depends on the 
angular variable 6 and a real variable x. 
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Inversion of the relation between the Wigner function and quadrature distributions becomes 
quite obvious, if we rewrite Eq. ( p.l5| ) in terms of the Fourier transforms of both the sides. The 
generating function for the quadrature distribution can be expressed using the Wigner function 
as: 

ZeiX) = /dxe^^Xx) 

= dx dye''^^W{xcos6 — ysm6,xsm6 + ycos6) 

= f dq fdpe'^^'''°'^-^P''''^^W{q,p). (3.17) 

The last expression is simply the Fourier transform of the Wigner function taken at the point 
(g cos 9, p sin 9). Thus, the projection relation expressed in terms of the Fourier transforms con- 
sists in the change of the coordinate system, from the Cartesian one (Wigner function) to the 
polar one (quadrature distributions). 

With this observation in hand, the way to invert Eq. ( [3.15D is straightforward: we need to write 



the Wigner function as the inverse Fourier transform, and to change the integration variables from 
Cartesian to polar. This allows us to insert the generating function for quadrature distributions: 

W{q,p) = ——r \\\d\ d^e-^''^™^^-*P^^'°^Ze(A). (3.18) 

(27r)^ y_oo Jo 

Expressing Zg(\) in terms of quadrature distributions and performing the integral over A yields: 

W{q,p) = ^ r dx r d9pe{x)^P \ — (3.19) 

where P denotes the principal value. The above formula is known as the inverse Radon trans- 
form. It is clearly seen that this transformation is singular. Therefore, its numerical implementa- 
tion is quite complicated. When processing experimental distributions pe{x), which are affected 
by statistical noise, one has to apply a regularization scheme. 

The close link between the quadrature distributions and the Wigner function could be noted 
already during the discussion of the balanced homodyne detector. The first expression for Zq{\) 
in Eq. ( |J!6| ) is exactly the symmetrically ordered characteristic function that appears in Eq. (p^, 
with C = iAe^Vv^. 

The first experimental realization of optical homodyne tomography has been demonstrated 
by Smithey et al. [0]. This seminal work has been followed by extensive theoretical and ex- 
perimental research. The effects of imperfect detection were analysed [^3]], and it was shown 
that in such a case the inverse Radon transform yields a generalized quasidistribution function 
with the ordering parameter equal to —(1 — z/)///, where r] is the efficiency of the photodetectors. 
Thus, detector losses result in blurring of the measured Wigner function. A more fundamental 
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problem related to optical homodyne tomography was the determination of other quantum state 
representations from homodyne statistics. In principle, once we have the Wigner function, we 
can evaluate the expectation value of any quantum observable O, and obtain for example the 
density matrix in the Fock basis. However, it would be appealing to reconstruct the observables 
directly from the homodyne statistics, in order to avoid the detour via the singular inverse Radon 
transform. The formula needed for this purpose is of the form: 

/oo rn 

dx / d9fo{x,9)pe{x), (3.20) 

oo Jo 

where /q(x, 9) is called the pattern function related to the observable O. The problem of deriving 
pattern functions for the elements of the density matrix in the Fock basis was studied first by 
D'Ariano et al. p5|]. It was later generalized to a more fundamental form [^ ^Sp. The statistical 



error of optical homodyne tomography has been thoroughly studied in a series of papers [|49|, |50 



51| , p20 . On the experimental side, optical homodyne tomography has been demonstrated for 
cw fields [|53|], and a beautiful gallery of squeezed states of light has been presented [[T0[]. An 
analogous tomographic method has been used to characterize transversal degrees of freedom of 
a laser beam [p4|]. 

3.4 Random phase homodyne detection 

In the context of optical homodyne tomography, a new technique for measuring light has been 
developed. This technique is balanced homodyne detection with the phase 9 made a uniformly 
distributed random variable [^. The distribution of events Pti{x) observed in such a case is 
described by phase-averaged homodyne statistics: 

Pn{x) = ^ f \9pe{x) = ^ T d9 {\x)ee{x\) . (3.21) 

In this regime, the phase sensitivity of homodyne detection is completely lost, and the phase- 
averaged homodyne statistics pn{x) contains information only on phase-independent properties 
of the measured light. Nevertheless, random phase homodyne detection has some advantages 
compared to direct photodetection. First, ultrafast sampling time can be achieved by using the 
local oscillator field in the form of a short pulse. Second, information on the photon distribution 
is carried by two rather intense fields, which can be detected with substantially higher efficiency 
than the signal field itself. This feature has enabled an experimental demonstration of even-odd 
oscillations in the photon distribution of the squeezed vacuum state [ ]56| ] . 

Let us now see, how the phase-averaged homodyne statistics depends on the photon distribu- 
tion. We shall use the fact that eigenvectors of the quadrature operator xe can be obtained from 
the position eigenvectors \x) by the unitary transformation \x)0 = e*^" "|x). This unitary trans- 
formation is diagonal in the Fock basis, and e*^" "|n) = e*"^|n). Introducing two decompositions 
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of the identity operator in the Fock basis, we have: 

Pnix) = 75- / ^6* ^ {\m){m\x)ee{x\n){n\) 

^^ -^0 m,n=0 

= y {m\x){x\n) / — e'("-'")^(|m)(n|) 

— ' /n Stt 

m,n=0 "^^ 

00 

= ^\{m\x)\^{\m){m\) (3.22) 

m=0 

Thus, occupations of the Fock states given by the expectation values (|m) (m|) contribute to the 
phase-averaged homodyne statistics with the coefficients 

\Hx)\' = ^^,Hl{x)e-\ (3.23) 

where Hm{x) denote Hermite polynomials. These coefficients correspond to the position distri- 
butions for the eigenstates of the harmonic oscillator. 

Integrating j97^(x) with appropriate pattern functions, we may reconstruct the photon statistics 
of the measured field, as well as other phase independent observables. An alternative method of 
processing the phase-averaged homodyne statistics, based on maximum-likelihood estimation, 
has been described in Ref . [ ^ . 
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Chapter 4 

Direct probing of quantum phase space 



Quasidistribution functions contain complete characterization of the quantum state. An inter- 
esting and nontrivial problem is how to determine quasidistributions from quantities which can 
be detected in a feasible experimental scheme. In the previous chapter, we have discussed two 
experimental techniques for measuring the quantum state of a light mode: double homodyne de- 
tection, and optical homodyne tomography. These techniques are based on detection of quadra- 
tures, which are continuous variables. The two-dimensional probability distribution observed in 
double homodyne detection yields directly the Q function of a light mode. In optical homodyne 
tomography, a family of one-dimensional projections of the Wigner function is measured and 
then processed numerically using the back-projection algorithm. 

In this chapter, we shall present a different approach to measuring quasidistributions of a 
light mode. We have seen in Chap. ^ that quasidistributions at a specific point of the phase 
space are given by expectation values of certain Hermitian operators. We shall demonstrate that 
this definition leads to a novel optical scheme for measuring quasidistribution functions of light. 
This scheme is based on photon counting. In contrast to the homodyne techniques discussed in 
the previous chapter, it is essential in our approach that the signal obtained from the detector is 
discrete, and that it is described by an integer variable characterizing the number of absorbed 
photons. 



Our starting point in Sec. ^J] will be a simple relation between the Wigner function and the 
photon statistics. We show that this relation can be implemented using a simple optical setup, 
which allows one to determine quasidistributions from photon statistics. In Sec. ^]2| we discuss 
the proposed scheme using the phase space picture. This picture explains in an intuitive way how 
physical parameters of the setup determine the measured quantity. In Sec. ^3| we generalize the 
relation linking the quasidistributions and the photon statistics. We also discuss effects of non- 
unit detector efficiency. The full multimode theory of the proposed setup is developed in Sec. B^. 



We show there that the direct method for measuring quasidistributions can be easily extended to 
multimode radiation. Finally, in Sec. Ol we discuss theoretically several examples of photon 



statistics which would be obtained from the photodetector when measuring quasidistribution 
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functions. 

4.1 Wigner function and photon statistics 

We will start from deriving a simple relation between the Wigner function at the origin of the 
complex phase space W{0) and the photon statistics. Let us take the operator ( [2.13[ ) for « = 
and s = 0, and expand it into a power series: 

W{0; 0) = - : exp (-2a^a) : 

TT 

n=0 



2 °° 
J2i-mn){n\. (4.1) 



TT 

n=0 



In the last step, we have used the normally ordered operator representation of the n photon 
number projection operator. The last expression contains a sum which assigns +1 to even Fock 
states, and —1 to odd Fock states. Therefore, the whole sum is simply the parity operator, and the 
Wigner function at the phase space origin is given, up to the front factor 2/n, by its expectation 



value [^. Taking the quantum average of Eq. ( |4~T| ) gives 



oo 



W{0) = {W{0; 0)) = - V(-l)>„ (4.2) 

TT -^^ 
n=0 

where the value p„ appearing in this expansion is just the probability of counting n photons by 
an ideal photodetector. Thus the photon statistics allows one to evaluate the Wigner function at 
the origin of the phase space. 

It would be appealing to generalize this relation to an arbitrary point of the phase space. In 
principle, the only thing we have to do is to shift the system or equivalently the frame of reference 
in the phase space. The problem is, how to realize this in practice for optical fields. We will show 
that this goal can be achieved using a very simple optical arrangement, and that in a certain limit 
the displacement transformation is realized. 



Let us consider the setup presented in Fig. ^^. We take the detected field to be a superposition 



of two single-mode fields, which we shall call the signal and the probe. The corresponding 
annihilation operators are denoted by as and dp respectively. The superposition is realized by 
means of a beam splitter BS with the power transmission characterized by the parameter T. 
In general, the action of the beam splitter is described by an SU(2) transformation between the 
annihilation operators of the incoming and outgoing modes p7|]. As the phase shifts appearing in 
this transformation can be eliminated by appropriate redefinition of the modes, the annihilation 
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Figure 4. 1 : Experimental setup for measuring directly quasidistribution functions of a single light mode. 
BS denotes the beam splitter, PD is the photodetector, and the annihilation operators of the modes are 
indicated. 



operator dout of the outgoing mode falling onto the detector surface can be assumed to be a 
combination 

dout = VTas - Vl - Tap. (4.3) 

The photon statistics of the field dout is used to evaluate the alternating series according to 
Eq. (|4.2D. In this way we obtain the Wigner function of the outgoing mode at the origin of the 
phase space. This quantity depends on the quantum state of both the modes as and dp. The 
Wigner function of the outgoing mode at the phase space origin is given in terms of the incoming 
modes by the expectation value of 

t(0; 0) = - : exp 



PVo, 



TT 

2 

IT 



^tXrint'^n 



exp 



-2r(d^ - y/{l-T)/Td^p){as - v/(l - T)/Tdp)) 



(4.4) 



This simple relation provides an interesting link between the detected quantity and the S mode. 
Let us consider the case when the probe field is a coherent state dp|a) = a\a) uncorrelated with 
the signal mode. Performing the quantum average over the P mode in Eq. ( PT^ is straightfor- 
ward due to the normal ordering of the operators. Taking the expectation value over the signal 
and recalling the normally ordered definition of quasidistributions ( [2.12[ ), we can recognize a 
quasidistribution of the signal field. Comparing its parameters, we obtain that the Wigner func- 
tion ( |4.2[ ) for the outgoing mode is proportional to an s = 1 — 1/T ordered quasidistribution 
function of the S mode: 



l^out(O) = -Ws 



'1-T 



T 



-a: 



1-T 



(4.5) 



taken at the point a/(1 — T)/Ta. Thus our setup delivers directly the value of the signal qua- 
sidistribution function at the phase space point dependent on the amplitude and the phase of the 
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probe coherent state. Since both these parameters can be controlled experimentally without dif- 
ficulties, we may simply scan the phase space by changing the amplitude and the phase of the 
probe field and thus determine the complete quasidistribution function. Eq. ( p3| ) shows that its 
ordering depends on the beam splitter transmission. For T near one the ordering is close to zero, 
which means that the measured quasidistribution approaches the Wigner function of the signal 
field. In contrast to optical homodyne tomography the Wigner function is measured directly and 
no sophisticated computer processing of the experimental data is necessary. The quantity mea- 
sured in the experiment is proportional to the quasiprobability distribution at the phase space 
point depending only on the amplitude and phase of the probe state. 

4.2 Phase space picture 

The quantity measured in the setup discussed in the previous section has an interesting phase 
space interpretation for arbitrary states of the S and P modes. To show this, we shall disentangle 
the two-mode operator defined in Eq. ( [4.4| ) using the following Gaussian integral of normally 
ordered operators for the S and P modes: 

W-out(0; 0) = ^ /"d'/5 : exp (-2(v^/3* - a^iVrp - ap)^ : 

X : exp (-2{VT^^(3* - al){VT^^P - as)) : . (4.6) 

Under the assumption that the S and P modes are uncorrelated, we can evaluate separately 
expectation values over the signal and the probe modes. It is easily seen that these expectation 
values yield the values of Wigner functions of the signal and the probe modes Ws{\/^ — Tf3) 
and Wp{\/Tl3). Thus we obtain the following expression for the quantity detected by our setup: 

W^out(O) = j d?l3Ws{Vl^Tl3)Wp{VTl3) 

= Y^j^'PWs{P)Wpiy/T/il - T)(3). (4.7) 

This formula establishes a connection between the photon number parity of the outgoing mode 
and the Wigner functions of the S and P modes. The object of interest in the above formula is 
the Wigner function of the signal field Ws{l3). The quantity we obtain from the measurement is 
VFout(O), which is given by the integral of Ws{(3) with the function Wp{^JT/{l -T)P). Let us 
now see, what information on the signal state can be obtained from this integral. 

In the case when the beam splitter splits the light equally, i.e. the power transmission is T = 
50%, we have a/T/(1 — T) = 1 and VFout(O) is simply a doubled overlap of the signal and probe 
Wigner functions. If we take the probe field to be a coherent state with variable amplitude and 
phase, we obtain the Q function of the signal field. Analogously to double homodyne detection 
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Wp{P) 




Figure 4.2: Phase space interpretation of the observable measured in the discussed setup. The Wigner 
function of the signal field Ws{/3) is integrated with the rescaled Wigner function of the probe field 

Wp{\/ Y^P). If the rescaling factor \ j^ is larger than one, the probe Wigner function becomes 
effectively contracted in all directions. Of course, the rescaled Wigner function does not describe physical 
fields at any point of the setup. 



discussed in Sec. p^ the probe Wigner function can be considered as a "quantum ruler", which 
smoothes the signal Wigner function to a positive definite phase space distribution. 

In a general case, the phase space parameterization of the probe Wigner function is rescaled 
by the factor ^JT/{l — T). It is easy to see that this factor can take an arbitrary positive value 
depending on the beam splitter transmission. The most interesting region is for T > 1/2, where 
the scaling factor ^JT /{I — T) is greater than one. This situation is shown pictorially in Fig. p!2 . 
The scaling factor effectively "contracts" the probe Wigner function in all the directions simul- 
taneously. Consequently, the area occupied by the rescaled probe Wigner function becomes 
smaller, and the integral in Eq. ( [4-.7[ ) provides more "local" information on the behaviour of the 
signal Wigner function Ws{l3). It is seen from the form of the scaling factor, that this effect of 
"contraction" grows unlimitedly with the beam splitter transmission tending to one. 

In a particular case of the coherent probe field \a), the P mode Wigner function is of the 
form 



Wp{(3) = -exp(-2|/5-ap). 



TT 



(4.8) 



When T tends to one, the rescaled probe Wigner function Wp{^JT/{l — T)P) approaches the 
form of the delta function, and the effect of smoothing becomes negligible. Thus, in the limiting 
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case of T ^ 1 the integral ( |4.7[ ) yields the value the signal Wigner function at a single phase 
space point. However, it it seen from Eq. ( ^3] ) that the rescaling of the probe Wigner function in 
the integral ( |4.7| ) has another consequence. When T tends to one, the factor multiplying the probe 
amplitude a becomes very small. This effect is seen also in Fig. p3| . Therefore, in order to scan 
the interesting region of the signal phase space we need to use a probe field of large intensity. 
From the mathematical point of view, the discussed limiting case involves two transitions: with 
T — i> 1 and |a| — > oo, such that the product a/(1 — T)/Ta is fixed. This product determines the 
point of the phase space at which the signal Wigner function is measured. 

Let us stress that the rescaled Wigner function, not obeying the Heisenberg uncertainty prin- 
ciple, does not describe any fields appearing physically in the setup; it is a purely abstract object 
introduced in the phase space interpretation of our measurement scheme. 

The integral representation derived in Eq. ( |4-.7D shows a connection of our setup with the 



model scheme of a phase space measurement, discussed in Sec. p^ . In this model scheme, a 
filter device — a "quantum ruler" — is introduced in addition to the system the measured phase 
space probability distribution is the convolution of the system and filter Wigner functions. Our 
scheme is more general, since the Wigner function of the filter can be rescaled by an arbitrary 
factor. Consequently the rescaled probe Wigner function does not have to obey the Heisenberg 
uncertainty principle and it may even approach the shape of a delta function, which leads to the 
direct measurement of the Wigner function. 

4.3 Generalization 

In the remaining discussion we will introduce two generalizations. First we will make our con- 
siderations more realistic by taking into account the imperfection of the photodetector. When the 
detector efficiency is i], the probability of counting n photons is given by the expectation value: 

p„ = ( : e-r^aLao J^^out^out)" \ _ ^^^^ 

The second extension is the substitution of the factor (—1)" inEq. ( [4-.2| ) by — (s + l)"/(s — 1)""^^, 
where s is a real parameter. The origin and the role of the parameters i] and s is different: r] 
describes experimental limitations, while s is an artificial number introduced in the numerical 
processing of the measured data. With these two parameters we obtain the following simple 
generalization of the formula ([4-.4D, when expressed in terms of the S and P modes 



2 ^ fs + iy _ „_„a,t„ae„, (^aLflout)" 






n=0 



2 / 2r/ 



vrfl — s) \ 1 — s 



exp I C^out'-^out 
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2r]T 



exp 



-{al - ^/{1-T)/Ta^p){as - ^{l^T)Jfap] 



7r(l — s) \ 1 — s 

(4.10) 

The third line of this equation suggests that the parameter s can be used to compensate the 
imperfectness of the photodetector. Indeed if we selected s = 1 — 77, we would determine the 
expectation value of : exp(— 2aoutaout) : regardless of the detector efficiency. But in this case the 
factor multiplying the probability of counting n photons is (1 — 2/77)" and its magnitude diverges 
to infinity with n ^ 00. Therefore we may expect problems with the convergence of the series. 
Even when the series is convergent, some singularities can be encountered in the processing 
of the experimentally measured photon statistics, which is affected by statistical fluctuations. 
Statistical noise may be a source of problems, as the increasing factor in Eq. ( }4-.10D causes that 
an important contribution comes from the "tail" of the experimental counts distribution, which 
usually has a very poor statistics. The diverging factor [(s + l)/(s — 1)]" amplifies fluctuations in 
this tail, and consequently the final result has a huge statistical error. The simplest way to avoid 
all these problems is to assume that the factors multiplying the counts statistics are bounded, 
which is equivalent to the condition s < 0. We shall discuss thoroughly the effects of statistical 
noise in the next chapter. This discussion will fully confirm the present conclusion drawn from 
qualitative arguments that in principle it is not possible to compensate for detector losses, and 
that in general we should restrict the range of s to nonpositive values. 

As before, let us consider the case when the coherent state \q) is employed as a probe. The 
expectation value of the generalized operator W^out (0; s) is again given by the quasidistribution 
function of the signal mode: 

(^iu?(0;s)) = —Ws I \l^^a--^——^ 1 . (4.11) 




Let us now analyze the ordering of this function. As we discussed in Sec. |Z4] , although from a 
theoretical point of view an arbitrarily ordered distribution contains the complete characterization 
of the quantum state, experimental errors make it difficult to compute higher ordered distributions 
from the measured one. Thus what is interesting is the highest ordering achievable in our scheme. 
Analysis of the role of the parameter s is the simplest, since the greater its value, the higher is 
the ordering obtained. Because we have restricted its range to nonpositive values, it should be 
consequently set to zero. Thus we are left with the product of the parameters 77 and T. For 
fixed 7] the highest ordering is achieved when T ^ 1, and its limit value is now —(1 — r^)///. 
Under the assumption that r] and T are close to one, the ordering of the measured distribution is 
effectively equal to this limiting value if the difference 1 — T is much smaller than 1 — r]. This 
is a realistic condition for currently existing photodetectors, which have the maximum efficiency 
about 80% [p9|]. Thus the highest ordering achievable in our scheme is effectively determined 
by the photodetector efficiency and is equal to —(1 — r/)///. It is noteworthy that this is exactly 
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equal to the ordering of the distribution reconstructed tomographically from data measured in the 
homodyne detection with imperfect detectors P3|]. 



4.4 Multimode approach 

So far, we have discussed the experimental scheme assuming that the signal and the probe beams 
are single-mode fields. We shall now present the full multimode theory of the direct scheme for 
measuring the quantum optical quasidistribution functions. In the single-mode description it was 
sufficient to use a pair of annihilation operators as and ap. The spatio-temporal characteristics 
of these two modes was not important in this approach, and it was implicitly assumed that the 
modes are matched perfectly at the beam splitter. We shall free our further analysis from these 
simplifying assumptions. 

Let us denote by 'EU (r, t) the positive-frequency part of the electric field operator at the 
surface of the detector. This field is a superposition of the signal and the probe fields combined 
at the beam splitter BS. Mathematical representation of this combination is a slightly delicate 
matter. If we wanted to express 'EU (r, t) in terms of the signal and probe field operators before 
the beam splitter, we would have to introduce appropriate propagators. This would obscure the 
physical picture of the measurement. Therefore we shall choose another notation for the signal 
and the probe fields, which will make the discussion much more transparent. We shall denote by 
Eg (r, t) the electric field operator of the signal beam that would fall onto the detector surface in 
the absence of the beam splitter BS. Analogously, let Ep (r, t) be the probe field at the detector 
surface, assuming that the beam splitter BS was replaced by a perfectly reflecting mirror. With 
these definitions, the field Eq^ (r, i) resulting from the interference of the signal and the probe 
beams is given simply by 

Ei+)(r,t) = VfEf\v,t) - v/r^TEi,+)(r,t). (4.12) 

We have assumed here that the characteristics of the beam splitter is constant over the spectral 
and polarization range of the considered fields. 

Further, we shall assume that the detected fields are quasi-monochromatic with the central 
frequency ujq. This will allow us to relate easily the number of photons to the energy of the 
field absorbed by the detector. Assuming that the direction of propagation of the field Eo^t (r, t) 
is perpendicular to the detector, the operator of the photon flux through the detector surface is 
given by 

iout = ^ / dt / d2rEl-V,^)EiutV,i) (4.13) 

"^0 J At J D 

where Ei^^ (r, t) = [Eo^t (r, t)]^^, and the temporal and the spatial integrals are performed respec- 
tively over the detector opening time At and its active surface D. The probability of registering 
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n photons is given by 



vJo J'n^-^y^ X ^ (4.14) 



Pn= { ■■ e 

n 

' s,p 



where r] is the detector quantum efficiency. 

As before, we will use the count statistics to calculate the average parity of the registered 
photons. It can be expressed in terms of the photon flux operator as: 

oo 

^(-1)>„ = ( :exp(-2r/Jo„0 : )s,p. (4.15) 

n=0 

If a coherent field is used as the probe, we may immediately evaluate the quantum expectation 
value over the P mode and obtain 

[VfE^s^\r,t) - v^r^TEp(r,t)]') :\ , (4.16) 

where Ep(r, t) = (Ep ^(r, t))p is the amplitude of the coherent probe field. 

We will now consider the signal field E^^ (r, t) in which a finite number of M modes is 
possibly excited. We shall denote the corresponding annihilation operators by Sj, and the mode 
functions by Uj(r, t), where i = 1,2,... M. Our goal will be to relate the photon statistics pn 
to the multimode quasidistribution characterizing these modes. Thus we decompose the signal 
field E5.+^(r,t) in the form 

M 



n=0 

X 



j=l 



where the operator V(r, t) is a sum of all the other modes remaining in the vacuum state. This 
part of the field does not contribute to the detector counts in the normally ordered expression 
given in Eq. ( [4-.16D , because its normally ordered moments are zero. 



Further, we shall assume that virtually all the excited part of the signal field is absorbed by 
the detector within the gate opening time. This allows us to write orthonormality relations for 
the mode functions Uj(r, t) in the form 

^ / dt / dhu*{r,t)uj{r,t) = 5ij (4.18) 

"-^0 J At J D 

where the integrals are restricted to the domain defined by the detection process. With these 
assumptions, we may simplify the exponent of Eq. ( |4.16p . It is convenient to introduce dimen- 
sionless amplitudes ctj, which are projections of the probe field onto the mode functions: 

2eoc f ^f ^o 



J M JD 



«. = -^ / dt / d^ru*(r,t)Ep(r,t). (4.19) 
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Using these amplitudes, we may write the measured quantity as: 

^(-1)>„ = ( : exp -2r/T^(al - ^(1 - T)/Tal){a - ^(1 - T)/Tc 



n=0 \ \ 1=1 




X exp 



-^"(^-^'(lfI*X^^^i^^<--«'i^-B-i'^) 



(4.20) 



The exponent appearing in the second line of the above expression results from the part of the 
probe field that is orthogonal (in the sense of Eq. ( 14.18D ) to the mode functions describing the 
excited component of the signal field. This exponent is equal to one if the probe field matches 
the M signal modes of interest. This condition can be written as: 



M 



Ep(r,t) = J]«,Ui(r,t). (4.21) 



i=l 



In this case, we can easily recognize in the quantum expectation value in Eq. ( [4.2(J| ) the multi- 
mode quasidistribution function defined in Eq. ( |2.24D , and write: 




E(-i»» H ^ M V^-. ■ ■ ■ . V^-'; -^ ■ («2) 



n=0 

Thus, the direct scheme allows one to measure multimode quasidistribution functions, even if 
the modes cannot be spatially separated. What one needs to do, is to combine the multimode 
signal field with appropriately chosen probe field, and to measure the count statistics of the 
resulting superposition. Let us note, that it is not necessary to resolve contributions to the count 
statistics from each of the modes; the only observable we need to reconstruct the multimode 
quasidistribution is the parity of the total number of photocounts. 

The direct scheme for measuring multimode quasidistributions can be also used if some of 
the modes are spatially separated. In this case, each mode (or group of spatially overlapping 
modes) has to be displaced in the phase space by combining at a beam splitter with a coherent 
probe field, and then measured using a photon counting detector. The parity of the number of 
photocounts obtained on all the detectors yields the value of the quasidistribution function at a 
point defined by the values of coherent displacements. 

4.5 Examples of photocount statistics 

We will close this chapter by presenting several examples of the photocount statistics for different 
quantum states of the signal probed by a coherent source of light. The most straightforward case 
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is when a coherent state \ao) enters through the signal port of the beam splitter. Then the statistics 
of the registered counts is given by the Poisson distribution: 

p\ao) ^ IHO^^-Jiao)^ (4_23) 

n! 



where J(ao) = ^^1/3— aoP is the average number of registered photons, and [3 = a/(1 — T)/Ta 
denotes the point of the phase space at which the value of the quasidistribution is measured. 
When the measurement is performed at the point where the quasidistribution of the signal field is 
centered, i.e., /? = ao, the fields interfere destructively and no photons are detected. In general, 
for an arbitrary phase space point, the average number of registered photons is proportional to 
the squared distance from ao- Averaging Eq. ( [4-.23[ ) over an appropriate P representation yields 



the photocount statistics for a thermal signal state characterized by an average photon number n: 

^" " (l + r7rn)"+i^" \n{l + 7]Tn)) ''''^ \l + vTn) ' ^"^'^"^^ 

where L„ denotes the nth Laguerre polynomial. 

A more interesting case is when the signal field is in a nonclassical state, which cannot de- 
scribed by a positive definite P function. Then the interference between the signal and the probe 
fields cannot be described within the classical theory of radiation. We will consider two nonclas- 
sical states: the one photon Fock state and the Schrodinger cat state. The most straightforward 
way to calculate the photocount statistics is to evaluate explicitly the general quasidistribution 
function, and to substitute as its parameters a/(1 — T)/Ta and — (1 — s — rjT)/rjT according to 
Eq. ( |4.11| ). Expanding this expression into the powers of (s + l)/(s — 1) yields the photocount 
statistics, which follows from the first line of Eq. ( |4.1(J| ). 

The photocount distribution for the one photon Fock state 1 1) can be written as an average of 
two terms with the weights rjT and 1 — rfT: 

pID = r^T[n - J[Q)f\dMlle-'^^) + (i _ ^T)M!e-''W . (4.25) 

TV. n\ 

The second term corresponds to the detection of the vacuum signal field. Its presence is a result 

of the detector imperfection and the leakage of the signal field through the unused output port 

of the beam splitter. This term vanishes in the limit of rjT -^ 1, where the Wigner function is 

measured in the setup. The first term describes the detection of the one photon Fock state. In 

Fig. ^j^(a) we show the statistics generated by this term for different values of /?. If the amplitude 

of the probe field is zero, we detect the undisturbed signal field and the statistics is nonzero only 

forn = 1. The distribution becomes flatter with increasing [3. Its characteristic feature is that it 

vanishes around n ^ J(0). 

For the Schrodinger cat state defined in Eq. ([2.18[) the photocount statistics is a sum of three 



terms: 






['^(^o)]\ -i(.o) + [-^(-«0)]\ -J(-»o) 



2(l + e-2|"oP) [ n\ n\ 
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(a) one photon Fock state 



(b) Schrodinger cat state 
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Figure 4.3: The photocount statistics of (a) the one photon Fock state and (b) the Schrodinger cat state for 
Of) = 3i, shown for several values of the rescaled probe field amplitude (3 = ^/{l — T)/Ta in the limit 
rjT = 1. 



CHAPTER 4. DIRECT PROBING OF QUANTUM PHASE SPACE 45 

. (4.26) 



_^2Re f [^^(/^* ~ "o)(/^ + O'o)]'' ^riT(a*0-aol3*)\ ^-{2-nT)\ao\^-riT\l3\^ 



The first two terais describe the two coherent components of the cat state, whereas the last one 
contributes to the quantum interference structure. In Fig. ^]^(b) we plot the photocount statistics 
for different values of P probing this structure, in the limit rjT ^ 1. The four values of P 
correspond to the cosine function in Eq. ( [2.19D equal to 1, 0, —1, and 0, respectively, for t = 0. It 
is seen that the form of the statistics changes very rapidly with /3. This behavior becomes clear 
if we recall that the Wigner function is given by the expectation value of the displaced photon 
number parity operator. Therefore, in order to obtain a large positive (negative) value of the 
Wigner function, the photocount statistics has to be concentrated in even (odd) values of n. 
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Chapter 5 
Practical aspects 



We will now discuss practical aspects of the direct scheme for measuring quasidistributions of 
a single light mode. An important problem is the role of the statistical noise. The proposed 
measurement scheme is based on the relation between the quasidistributions and the photocount 
statistics. In a real experiment the statistics of the detector counts cannot be known with perfect 
accuracy, as it is determined from a finite sample of N measurements. This statistical uncertainty 
affects the experimental value of the quasidistribution. Theoretical analysis of the statistical error 
is important for two reasons. First, we need an estimate for the number of the measurements 
required to determine the quasidistribution with a given accuracy. Such an estimate is needed 
when designing an experiment. The total number of measurements is usually limited by various 
factors, such as temporal stability of the optical setup. Estimation of the statistical error tells us. 



how precise result can be expected in a realistic scheme. Secondly, we have seen in Sec. ^3| that 
one may attempt to compensate the imperfection of the detector and the non-unit transmission 
of the beam splitter by appropriate numerical processing of the measured statistics. However, a 
preliminary qualitative discussion of the statistical noise strongly indicated that such a procedure 
may amplify the statistical noise. Our present calculations will provide a detailed, quantitative 
analysis of this problem. In Sec. |5TT| we define the observable whose statistical properties will be 
studied. In Sec. ^]2| we derive expressions for the mean value and the statistical variance. The 
basic element of this derivation is the multinomial distribution, which defines the probability of 
obtaining a specific photocount statistics from a series of N measurements. The behaviour of 
the mean value and statistical variance is discussed in Sec. ^]3|. Particular attention is paid to 
the possibility of loss compensation. We present a collection of pathological cases, where the 
compensation leads to an explosion of the statistical uncertainty. These cases show clearly that 
in a general case the compensation of detector losses is not possible. 

In Sec. ^]^ we discuss effects of the mode-mismatch between the signal and the probe fields. 
In practice, the fields superposed at a beam splitter never exhibit 100% visibility of interfer- 
ence. This fact has to be taken into account, when we analyse the result of a real experiment. 
The effects of the mode-mismatch will be discussed using the multimode theory developed in 

47 
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Sec. HI. 



5.1 Photon count generating function 

In order to make our discussion more transparent, we shall redefine the quantity evaluated from 
the photocount statistics p„ to the form: 

n=0 ^ ^ 

Thus, n(s) corresponds to the observable defined in Eq. ( |4.1(J| ) without the front normalization 
factor. We shall call n(s) the photon count generating function (PCGF), as the full photon 
statistics can be retrieved from n(s) as an analytical function of s. In our scheme, the PCGF is 
given by the expectation value 




n(.) = (:expl-^^l :), (5.2) 

where jTout is the operator of the time-integrated flux of the light incident onto the surface of the 
detector. This operator can be expressed in terms of the signal and probe fields as 

Jout = {VTal - Vl-Ta^p){VTas - Vl-Tap), (5.3) 

with T being the beam splitter power transmission. When a coherent state \a)p is used as a 
probe, we have 

and the PCGF is proportional to the quasidistribution function of the signal mode. 

5.2 Statistical error 

In a real experiment, the photon statistics is obtained from a finite series of A^ measurements. The 
result of these measurements has the form of a histogram {fc„}, where fc„ denotes the number 
of measurements when n photons have been detected. Dividing kn by the total number of mea- 
surements, we obtain an estimate for the photon distribution p„. This estimate is subsequently 
used to evaluate the PCGF defined in Eq. (p.lj). In the analysis of the statistical properties of the 
PCGF, we shall introduce a cut-off parameter K for the maximum photon number. In this way, 
the analysed quantity will depend on a finite number of variables /c„, where n = 0,1, . . . , K. 
Thus, we consider an experimental estimate for 11 (s) in the form: 

n=0 ^ ^ 
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Our goal is to see how well this estimate approximates the ideal quantity n(s). Due to 
extreme simplicity of the relation between the count statistics the quasidistributions, it is possible 
to perform a rigorous analysis of the statistical error and to obtain an exact expression for the 
uncertainty of the final result. The basic tool in our analysis is the probability distribution for 
the histograms {A;„} which can be obtained from N experimental runs. Because all the runs are 



statistically independent, the set of fc„s obeys the multinomial distribution [ pO| ] 



V{ko,ki,...,kK) 



N\ 



K 



^-En=ofcn 



K 



kom...(N-Zn=okn 



-P0'"P^...P^- l-$^Pn 



n=0 



(5.6) 

In this formula, p„ is the ideal, noise-free photocount distribution that would be obtained in the 
limit of the infinite number of measurements. 

Using the distribution given in Eq. (|5.6D, we may calculate quantities characterizing statistical 
properties of the experimental PCGF. In order to see how well nexp(s) approximates the ideal 
quantity we will find its mean value and its variance. This task is quite easy, since the only 
expressions we need in the calculations are the following moments: 



Kr. 



klkji 



Npn, 

N{N - l)pipn + 6lnNpn. 



We use the bar to denote the statistical average with respect to the distribution V{kQ, 
Given this result, it is straightforward to obtain: 
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(5.7) 
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(5.8) 



and 
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(5.9) 



The error introduced by the cut-off of the photocount statistics can be estimated by 



|nexp(s) - n(s) 



n=K+l 



S + 1 
S- 1 
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n=K+l 
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(5.10) 



The variance ^Ilg^p, derived in Eq. (|5^), is a difference of two terms. The second one is simply 
the squared average of Ilexp. The first term is a sum over the count statistics multiplied by the 
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powers of a positive factor [{s + l)(s — 1)]^. If s > 0, this factor is greater than one and the sum 
may be arbitrarily large. In the case when the contribution from the cut tail of the statistics is 
negligible, i.e., if fC ^ oo, it can be estimated by the average number of registered photons: 

E 731 Pn>l + TJ-^{vJo.). (5.11) 

n=0 ^ ^ ^ ^ 

Thus, the variance grows unlimited as we probe phase space points far from the area where the 
quasidistribution is localized. Several examples in the next section will demonstrate that the 
variance usually explodes much more rapidly, exponentially rather than linearly. This makes 
the compensation of the detector inefficiency a very subtle matter. It can be successful only for 
very restricted regions of the phase space, where the count statistics is concentrated for a small 
number of counts and vanishes sufficiently quickly for larger n's. 

Therefore, in order to ensure that the statistical error remains bounded over the whole phase 
space, we have to impose the condition s < 0. Since we are interested in achieving the highest 
possible ordering of the measured quasidistribution, we should consequently set s = 0. For 
this particular value the estimations for the uncertainty of Ilexp take a much simpler form. The 
error caused by the cut-off of the count distribution can be estimated by the "lacking" part of the 
probability: 

K 

|nexp(o)-n(o)|<i-5^p„, (5.12) 

n=0 

which shows that the cut-off is unimportant as long as the probability of registering more than K 
photons is negligible. The variance of Ilexp is given by 



^ne%(o) = ^ 



■ K 

5^Pn-(nexp(0) 

n=0 



1 
< — 



nexp(o) 



< ^. (5.13) 



Thus, the statistical uncertainty of the measured quasidistribution can be simply estimated as 
1/Vn multiplied by the proportionality constant given in Eq. (^^. It is also seen that the 
uncertainty is smaller for the phase space points where the magnitude of the quasidistribution is 
large. 

5.3 Compensation of detector losses 

We will now consider several examples of the reconstruction of the quasidistributions from the 
data collected in a photon counting experiment. Our discussion will be based on Monte Carlo 
simulations compared with the analytical results obtained in the previous section. 

First, let us note that the huge statistical error is not the only problem in compensating the 
detector inefficiency. If s > 0, the sum ( |5.8| ) does not even have to converge in the limit of 
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K ^ oo. An example of this pathological behaviour is provided by a thermal state, which has 
been calculated in Eq. ( [4-.24[ ). For the zero probe field we obtain 

which shows that if s > 0, then for a sufficiently intense thermal state the magnitude of the 
summand is larger than one and consequently the sum diverges, when K ^ oo. This behaviour 
is due to the very slowly vanishing count distribution for large n, and it does not appear for the 
other examples of the count statistics derived in Sec. ^3[ 

In Fig. |5Tl] we plot the reconstructed quasidistributions for the coherent state \ao = 1) and 
the one photon Fock state. Due to the symmetry of these states, it is sufficient to discuss the 
behaviour of the reconstructed quasidistribution on the real axis of the phase space. The cut-off 
parameter is set high enough to make the contribution from the cut tail of the statistics negligibly 
small. The quasidistributions are determined at each phase space point from the Monte Carlo 
simulations of A^ = 1000 events. The grey areas denote the statistical uncertainty calculated 
according to Eq. (|5^). The two top graphs show the reconstruction of the Wigner function in 
the ideal case r]T = 1. It is seen that the statistical error is smaller, where the magnitude of the 
Wigner function is large. In the outer regions it approaches its maximum value 1/^/N. The effect 
of the nonunit r]T is shown in the center graphs. The measured quasidistributions become wider 
and the negative dip in the case of the Fock state is shallower. In the bottom graphs we depict 
the result of compensating the nonunit value of r]T by setting s = 1 — r]T. The compensation 
works quite well in the central region, where the average number of detected photons is small, 
but outside of this region the statistical error explodes exponentially. Of course, the statistical 
error can be in principle suppressed by increasing the number of measurements. However, this 
is not a practical method, since the statistical error decreases with the size of the sample only as 

1/Vn. 

The reconstruction of the interference structure of the Schrodinger cat state is plotted in 
Fig. P?^ . We have used the state defined in Eq. ( [2.18D with ao = Si. The interference structure is 



very fragile, and its precise measurement requires a large sample of events. In the case of the pre- 
sented plot, N = 5000 simulations were performed at each phase space point. Comparison of the 
top and the center graphs shows how even relatively small imperfection destroys the interference 
pattern. The data collected in a non-ideal setup can be processed to recover the Wigner function, 
but at the cost of a significantly larger statistical error, as it is shown in the bottom graph. Outside 
the interference structure, we again observe the exponential explosion of the dispersion due to 
the increasing intensity of the detected light. 

From the above examples we clearly see, that compensation of any efficiency lower than 
100% is a delicate matter. This conclusion may be quite puzzling, if we compare it with recently 
published results on compensation of losses in photodetection [[T3p. It was shown there that 
the true photon distribution can be reconstructed from data measured by an imperfect detector. 
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Figure 5.1: Reconstruction of the quasiprobability distributions of the coherent state |ao = 1) (left) and the 
one photon Fock state (right) from N = 1000 events. The solid lines are the analytical quasidistributions 
and the grey areas mark the statistical dispersion. The plots are parameterized with the rescaled probe 
field amphtude /3 = y^(1 - T)/Ta. 
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Figure 5.2: Reconstruction of the interference structure of the Schrodinger cat state for uq = 3i from 
N = 5000 events at each point. 
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Figure 5.3: Reconstruction of the Wigner function of the coherent state and the one photon Fock state 
from the count statistics cut at K = 11, for r]T = 1 and s = 0. The number of events is A^ = 10^. 



provided that its efficiency is greater than 50%. If this condition is fulfilled, the photon distribu- 
tion can be obtained from experimental data via the so-called inverse Bernoulli transformation, 
and the statistical properties of such reconstruction have been shown to behave regularly. One 
could think of using this recipe in our scheme for measuring quasidistribution functions: first, 
application of the inverse Bernoulli transformation would yield the loss-free photon distribution, 
and then evaluation of the alternating series would yield the unblurred Wigner function. How- 
ever, a simple calculation shows that such a two-step method is completely equivalent to setting 
s = 1 — r]T in the PCGF, and we end up with all the difficulties discussed above. Explanation of 
this seeming contradiction is simple: quantities obtained from the inverse Bernoulli transforma- 
tion have strongly correlated statistical errors, which accumulate when evaluating the alternating 
series. We shall discuss this issue in detail in Sec. [7^ , using the general theory of statistical 
uncertainty in photodetection measurements. 

Finally, let us look at the effect of cutting the statistics at a finite value. Fig. ^3| shows 
the Wigner functions for the one photon coherent and Fock states reconstructed from the count 
distributions cut at K = 11. We performed a large number of A^ = 10^ simulations in order to 
get the pure effect of the cut-off that is not spoiled by the statistical uncertainty. The grey areas 
show the cut-off error, estimated using Eq. ( |5.12D . The reconstruction works well as long as the 
probability of detecting more than K photons is negligible. When the average number of incident 
photons starts to be comparable with the cut-off parameter, "ghost" structures appear. When we 
go even further, the Wigner function again tends to zero, but this is a purely artificial effect due 
to the virtually vanishing count distribution below K. This kind of "ghost" structures resulting 
from the cut-off were observed in the characterization of the motional state of a trapped ion [[TT|]. 
In this experiment, the analog of photon statistics Pn was reconstructed in an indirect way by 
monitoring the fluorescence [^, and it was necessary to truncate the statistics p„ in order to 
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keep the error on a reasonably low level. We have seen that in the regions of the phase space 
where a substantial part of the statistics pn is lost by the truncation, we do not have sufficient 
data to determine the value of the Wigner function. 

5.4 Mode mismatch 

Another experimental imperfection that occurs in a realistic setup is the non-unit matching of the 
modes interfered at the high-transmission beam splitter BS. In order to analyse consequences of 
the mode mismatch, we will use the multimode approach developed in Sec. ^i^. Let us take the 
signal field and the coherent probe field to be of the form 

Ep(r,t) = aup{r,t) (5.15) 

where us{r, t) and up(r, t) are the corresponding normalized mode functions, and the operator 
V(r, t) is the sum of all other signal modes remaining in the vacuum state. We will assume that 
the support of these mode functions lies within the domain defined by the detector surface and 
the gate opening time. This allows us to write the normalization of the mode functions as 

dt /d2r|u5(r,t)|2= / dt /d2r|up(r,t)p = |^ (5.16) 

At JD JAt Jd ^^oC 

Further, we will assume that the overlap of the functions U5'(r, t) and up(r, t) is real and positive. 
This can be always achieved by multiplying up(r, t) and a by appropriate conjugated phase 
factors. We will denote 

^ [ dt [ d\u*sir,t)up{T,t) = Ve. (5.17) 

Under these assumptions, we may represent the PCGF calculated from the count statistics 
given by Eq. (|4.14]) in the form: 



n(s) = /: exp (-^^[Talas - ^/^T{1 - T){ala + asa*) + (1 - T)\am :\ (5.18) 
Rearranging the terms in the exponent yields: 



^(1-^)t.. f./^(l-^)„,. ^~'~vT\_f 2r/(l-r)(l-0,„,,2 



(5.19) 

Thus, if the signal and the probe modes are not matched perfectly, the PCGF is given by the 
quasidistribution function of the signal, but multiplied by a Gaussian envelope exp[— 27/(1 — 
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T)(l — ^)\a\'^ / (1 — s)]. This envelope is centered at the origin of the phase space, and the faster 
it decays, the larger is the mode mismatch, characterized by the difference 1 — ^. Additionally, 
the parameterization of the signal quasidistribution is rescaled by v^. 

The effect of the mode-mismatch is more severe in outer regions of the phase space. It is 
not important if the quasidistribution is localized around the center of the phase space within the 
width of the Gaussian envelope. In particular, for the vacuum signal state we obtain: 

U{s) = exp f-^n^l^p-laA . (5.20) 

This result does not depend at all on the mode overlap parameter ^, simply because with the vac- 
uum signal field no interference occurs at the beam splitter BS and all the recorded photons come 
from the probe beam. The above expression can be interpreted as the quasidistribution function 
of the signal mode characterized by the mode function up(r, t), which perfectly overlaps with 
the probe field. 

In a general case, we may rewrite Eq. ( |5.19[ ) to the form: 



liis) = •': S' Ws\ l/ ^'^^^' a; 





l^iliL^a,— ^1. (5.21) 



where W^a^dP', s) is the vacuum quasidistribution function. This representation allows us to in- 
terpret Eq. ( |5.19| ) as a two-mode quasidistribution which is a product of the signal mode qua- 
sidistribution and an additional vacuum quasidistribution. These two modes are probed at the 
phase space points proportional to ^/^a and ^/T^^a. The amplitude y/^a describes the part 
of the probe field that overlaps perfectly with the signal field, whereas ^/T^^a corresponds to 
the orthogonal remainder (orthogonality is understood here in the sense of the scalar product 
between the mode functions). 

The parameter ^ characterizing the overlap of the signal and the probe modes can be related to 
the visibility of the interference. This expression will be useful in the discussion of the practical 
realization of the scheme. Let as assume that the signal mode is in a coherent state |ao)- The 
intensity of the field measured by the photodetector is given by 



/ = |yr^Ta- v^ao|' + (l-0^|ao|', (5.22) 

where a is the amplitude of the coherent probe mode. It is seen that for given ao the minimum 



intensity is obtained for y/1 — Ta = \/^ao, and it equals to: 

/min = (1 - OT\ao\^. (5.23) 
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When we now vary with the phase of a, keeping its absolute value fixed, the maximum intensity 
is achieved for Vl — Ta = —-/^aQ, and its value is: 

J^ax = 4eT|ao|' + (1 - i)T\a^?. (5.24) 

A simple calculation shows that the interference visibility v can be expressed using the overlap 
parameter,^ as: 

-'max -'min ^S zr- oc-\ 

V = = . (5.25) 

I 4- I ■ 2-\-2f 

-'max T^ -'mm -^ T^ -^S 

Inverting this relation, we obtain that ^ = v/{2 — v). 

Finally, let us note that the effect of mode-mismatch in our scheme is quite different from 
balanced homodyne detection, where it can be simply included in the overall detection efficiency 
parameter. In our scheme, it generates a Gaussian envelope multiplying the measured quasidis- 
tribution function, and its importance depends on the probed point of the phase space. 
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Chapter 6 
Experiment 



We shall now present direct measurement of the Wigner function by photon counting, using the 
method presented in Chap. ^. Previous measurements of the Wigner function of light, performed 
at the University of Oregon [^ and Universitat Konstanz [|T^, were realizations of optical ho- 
modyne tomography. In our experiment, the Wigner function is determined directly from the 
statistics of photocounts. Apart from the different principle of the measurement, we use a dif- 
ferent technique for light detection. As we discussed in Chap. ^ optical homodyne tomography 
is based on detection of the signal light superposed on a strong, classical local oscillator. The 
light incident on photodetectors has macroscopic intensity, and it is converted into an electronic 
current with the help of p-i-n photodiodes. At this level of intensity, it is not possible to resolve 
contributions from single photons, and the electric current is practically a continuous variable. 
Information on the quantum state of the measured light is contained in fluctuations of the differ- 
ence signal between two detectors. This signal is recorded using an analog-to-digital converter. 

In our scheme, the intensity of the detected light is comparable with the intensity of the signal 
field itself. Therefore, we need to use a detector that is sensitive to single optical photons. Cur- 
rently, the most efficient commercially available detectors on single-photon level are avalanche 
photodiodes operated in the so-called Geiger mode. The Geiger mode of operation consists in 
biasing the diode slightly above breakdown. Absorption of a photon triggers the breakdown, 
which is a macroscopic, recordable event. The diode is placed in a circuit which quenches the 
breakdown by lowering the voltage, and after a while restores the higher bias thus preparing 
the diode for the detection of a next photon. The electronic signal obtained from the diode has 
the form of pulses of uniform height and duration, which correspond to single detection events. 
These pulses, after shaping, can be counted using a standard digital logic device. 

We open this chapter with a review of the principle of the measurement in Sec. ^71]. The 
experimental setup is described in Sec. E^ and the results of the measurements are reported in 



Sec. |0| . In Sec. ^^^ we discuss the effect of various experimental factors, and summarize the 



presentation of the experiment. 
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6.1 Principle 

In order to make the presentation of the experiment self-contained, let us start with a brief review 
of the principle of the measurement. The Wigner function at a given phase space point is itself a 
well defined quantum observable. Furthermore, the measurement of this observable can be im- 
plemented for optical fields using an arrangement employing an auxiliary coherent probe beam. 
The amplitude and the phase of the probe field define the point in the phase space at which the 
Wigner function is measured. This allows one to scan the phase space point-by-point, simply by 
changing the parameters of the probe field. 

Our experiment is based on the representation of the Wigner function at a complex phase 
space point denoted by a as the expectation value of the following operator: 

2 °° 
^(«) = - V(-l)"^(a)|n)(n|Dt(a), (6.1) 

n ^ — ^ 

n=0 

where D(a) is the displacement operator and \n) denote Fock states, h\n) = n\n). Thus, W{a) 
has two eigenvalues: 2/n and — 2/7r, corresponding to degenerate subspaces spanned respec- 
tively by even and odd displaced Fock states. Practical means to translate this formula into an 
optical arrangement are quite simple. The displacement transformation can be realized by super- 
posing the measured field at a low-reflection beam splitter with a strong coherent probe beam. 
The value of the displacement a is equal in this setup to the reflected amplitude of the probe 
field. Furthermore, the projections on Fock states can be obtained by photon counting assuming 
unit quantum efficiency. These two procedures, combined together, provide a practical way to 
measure the Wigner function at an arbitrarily selected phase space point a. 

6.2 Setup 

The experimental setup we used to measure the Wigner function is shown schematically in 
Fig. p] . In principle, it is a Mach-Zender interferometric scheme with the beams in two arms 
of the interferometer serving as the signal and the probe fields. An attenuated, linearly polarized 



(in the plane of Fig. ^Jj) 632.8 nm beam from a frequency- stabilized single-mode He:Ne laser 
is divided by a low-reflection beam splitter BSl. The weak reflected beam is used to generate 
the signal field whose Wigner function will be measured. The state preparation stage consists 
of a neutral density filter ND and a mirror mounted on a piezoelectric translator PZT. With this 
arrangement, we are able to create pure coherent states with variable phase as well as their in- 
coherent mixtures. Though these states do not exhibit nonclassical properties, they constitute a 
nontrivial family to demonstrate the principle of the method, which provides complete charac- 
terization of both quantum and classical field fluctuations. 

The strong beam leaving the beam splitter BSl plays the role of the probe field with which 
we perform the displacement transformation D{a). In order to scan the phase space one should 
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Figure 6.1: The experimental setup for measuring the Wigner function. BSl and BS2 are quartz plates 
serving as low-reflection beam splitters. The quantum state is prepared using the neutral density filter ND 
and a mirror mounted an a piezoelectric translator PZT. The electrooptic modulators EOMl and E0M2 
control respectively the amplitude and the phase of the point at which the Wigner function is measured. 
The signal field, after removing spurious reflections using the aperture A, is focused on a single photon 
counting module SPCM. 



be able to set freely its amplitude and phase, which define respectively the radial and angular 
coordinate in the phase space. The amplitude modulation is achieved with a half-wave plate, 
a longitudinal Pockels cell EOMl, and a polarizer oriented parallel to the initial direction of 
polarization. The phase modulation is done with the help of an ADP crystal electrooptic phase 
modulator E0M2 on the signal field. This is completely equivalent to modulating the probe field 
phase, but more convenient for technical reasons: in this arrangement optical paths in both the 
arms of the Mach-Zender interferometer are approximately the same, and better overlap of the 
signal and the probe modes is achieved at the output of the interferometer. 

The signal and the probe fields are interfered at a nearly completely transmitting beam splitter 
BS2 with the power transmission T = 98.6%. In this regime, the transmitted signal field effec- 
tively undergoes the required displacement transformation. Spurious reflections that accompany 
the beam leaving the interferometer are removed using the aperture A. Finally, the transmitted 
signal is focused on an EG&G photon counting module SPCM-AQ-CD2749, whose photosen- 
sitive element is a silicon avalanche diode operated in the Geiger regime. The overall quantum 
efficiency of the module specified by the manufacturer is r^ > 70%. The count rate is kept low 
in the experiment, so that the chance of two or more photons triggering a single avalanche signal 
is very small, and the probability of another photon arriving during the detector dead time can 
be neglected. Under these assumptions, each pulse generated by the module corresponds to the 
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detection of a single photonQ The pulses are acquired by a computer, which also controls the 
voltages applied to the electrooptic modulators. The interference visibility in our setup has been 
measured to be f > 98.5%, and the phase difference between the two arms was stable up to few 
percent over times of the order of ten minutes. 

6.3 Results 

The voltages applied to electrooptic modulators were generated by high- voltage power supplies 
controlled by analog output ports of a multifunction I/O card (National Instruments PCTMIO- 
16E-4). A typical scan of the phase space consisted of sampling a sequence of circles with 
increasing radius. This was because changing the voltage applied to the amplitude modulator 
required a settling time of the order of 1 s. For a fixed amplitude the phase could be scanned 
much faster, as the phase modulator driver had the bandwidth up to 10 kHz. For a selected point 
of the phase space, the photon statistics was collected using a digital counter on the same I/O 
card, operated in the buffered event counting mode. 

In Fig. |6|^ we depict the measured Wigner functions of the vacuum, a weak coherent state, 
and a phase diffused coherent state. Phase fluctuations were obtained by applying a 400 Hz sine 
waveform to the piezoelectric translator. For all the plots, the phase space was scanned on a grid 
defined by 20 amplitudes and 40 phases. The scaling of the radial coordinate is obtained from 
the average number of photons rivac detected for the blocked signal path. Thus the graphs are 

1/2 

parameterized with the complex variable f3 = e^'^nj^c , where if is the phase shift generated by 
the phase modulator E0M2. At each selected point of the phase space, the photocount statistics 
Pn{P) was determined from a sequence of A^ = 8000 counting intervals, each r = 30/US long. 
The duration of the counting interval r defines the temporal envelope of the measured mode. The 
count statistics was used to evaluate the alternating sumQ 

oo 

n(/5) = $^(-i)>n(/?), (6.2) 

n=0 

which, up to the normalization factor 2/ii is equal to the Wigner function of the measured state. 
Statistical variance of this result can be estimated by Var[n(/3)] = {1 — [n(/5)]^}/A^, according 
to the discussion in Sec. ^^. Thus, the statistical error of our measurement reaches its maximum 
value, equal to l/N^^'^ ^ 1.1%, when the value of the Wigner function is close to zero. 



'The operating mode of the SPCM resuhs in a certain amount of extraneous pulses originating from dark counts 
and afterpulsing. The dark count rate of our module is less than 100 s^'^. This gives on average < 3 • lO^"^ during a 
single counting interval, which is 30 /is long. The typical afterpulsing probability is 0.2%. Thus, both these effects 
give a negligible contribution to the measured count statistics. 

^Here we explicitly write dependence of 11 on the phase space point /3. The number s used in the previous 
chapter as a parameter of 11 is now fixed and equal to zero. The operational phase space point (3 can be expressed 
by the amplitude a of the probe field as (3 — ^ri{l — T)a. 
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Figure 6.2: The measured Wigner functions of (a) the vacuum, (b) a weak coherent state, and (c) a phase 
diffused coherent state. The contour plots depict interpolated heights given by mutliples of 0. 1 for the 
plots (a) and (b), and by 0.08, 0.14, 0.20, 0.26, 0.32 for the plot (c). 

The Wigner functions of the vacuum and of the coherent state are Gaussians centered at the 
average complex amplitude of the field, and their widths characterize quantum fluctuations. It 
can be noticed that the measured Wigner function of the coherent state is slightly lower than 
that of the vacuum state. In the following, when discussing experimental imperfections, we shall 
explain this as a result of non-unit interference visibility. In the plot of the Wigner function of 
the phase diffused coherent state, one can clearly distinguish two outer peaks corresponding to 
the turning points of the harmonically modulated phase. Another set of results is presented in a 
colour plate at the end of this chapter. 

6.4 Discussion 



There are several experimental factors whose impact on the result of the measurement needs 
to be analyzed. First, there are losses of the signal field resulting from two main sources: the 
reflection from the beam splitter BS2 and, what is more important, imperfect photodetection 
characterized by the quantum efficiency i]. Analysis of these losses performed in Sec. ^3] shows, 
that in such a case the alternating series evaluated from photocount statistics is proportional to 
a generalized, s-ordered quasidistribution function W{a; s), with the ordering parameter equal 
s = —(1 — rjT) /rjT. 



64 MEASURING QUANTUM STATE IN PHASE SPACE 

In addition, the two modes interfered at the beam splitter BS2 are never matched perfectly. 
The effects of the mode mismatch have been analysed in Sec. ^j^. We will now apply this 
analysis to the experimental results, but first let us recall briefly the physical picture of mode 
mismatch. For this purpose, we need to consider the normalized mode functions describing the 
transmitted signal field and the reflected probe field. The squared overlap ^ of these two mode 
functions can be related to the interference visibility f as ^ = f /(2 — t>). In order to describe the 
effects of the mode mismatch, we need to decompose the probe mode function into a part that 
precisely overlaps with the signal, and the orthogonal remainder. The amplitude of the probe 
field effectively interfering with the signal is thus multiplied by ^^/^, and the remaining part of 
the probe field contributes to independent Poissonian counts with the average number of detected 
photons equal (1 — ^) |/3p. Consequently, the full count statistics is given by a convolution of the 
statistics generated by the interfering fields, and the Poissonian statistics of mismatched photons. 
A simple calculation shows, that the alternating sum evaluated from such a convolution can be 
represented as a product of the contributions corresponding to the two components of the probe 
field: 

n(/?) = exp[-2(l-OI/3r] 

Here on the right-hand side we have made use of the theoretical results for imperfect detection 
obtained in Sees. ^3] and |S!4] . Specializing the above result to a coherent signal state |ao) with 
the amplitude ao, yields: 



n(/3) = exp[-2|/5 - V^«o|' - 2(1 - OvT\ao\'']. (6.4) 



Thus, in a realistic case n(/5) represents a Gaussian centered at the attenuated amplitude y/^rjTao, 
and the width remains unchanged. This Gaussian function in multiplied by the constant factor 
exp[— 2(1 — ^)r]T\ao\'^]. For our measurement, ^ Ri 97% and r]T\ao\'^ ^ 1.34, which gives the 
value of this factor equal 0.92. This result agrees with the height of the experimentally measured 
Wigner function of a coherent state. 

In Fig. |5|3| we compare the phase-averaged measured Wigner function for the vacuum state 
with theoretical predictions. Agreement between the experimental points and the Gaussian curve 
is very good. This plot can be used to estimate the amount of excess thermal noise in the laser 
radiation. Let us assume that from the average nyac = |/3p registered photons a constant fraction 
Xth originates from thermal noise. The P-representation of the field that is effectively detected is 
given by 

^n' = r^ exp -— . (6.5) 
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Figure 6.3: Comparison of the phase-averaged experimental Wigner function for the vacuum state (points 
with error bars) with the theoretical prediction (line), given by the Gaussian exp(— 2|/J|^). 



The average photon number parity measured for such a field reads: 

1 / 2(1 



U{i3) 



exp 



Xth 



(6.6) 



2x.h|/3P + l ^V 2x,h|/3P + l 
and in the presence of thermal noise exhibits departure from a pure Gaussian shape. The perfect 



agreement observed in Fig. Q3]confirms that the contribution of thermal noise is negligibly small. 

Concluding, let us compare the demonstrated direct method for measuring the Wigner func- 
tion with the optical homodyne tomography approach. An important parameter in experimental 
quantum state reconstruction is the detection efficiency. Here better figures are exhibited by the 
homodyne technique, which detects quantum fluctuations as a difference between two rather in- 
tense fields. Such fields can be efficiently converted into photocurrent signals with the help of 
p-i-n diodes. It should be also noted that an avalanche photodiode is not capable of resolving 
the number of simultaneously absorbed photons, and that it delivers a signal proportional to the 
light intensity only in the regime used in our experiment. However, continuous progress in sin- 
gle photon detection technology gives hope to overcome current limitations of photon counting 
[p9|]. Alternatively, the displacement transformation implemented in the photon counting tech- 
nique can be combined with efficient random phase homodyne detection. This yields the recently 
proposed scheme for cascaded homodyning [|62|]. 

The simplicity of the relation ( |6.1[ ) linking the count statistics with quasidistribution func- 
tions allows one to determine the Wigner function at a given point from a relatively small sample 
of experimental data. This feature becomes particularly advantageous, when we consider detec- 
tion of multimode light. Optical homodyne tomography requires substantial numerical effort to 
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reconstruct the multimode Wigner function. In contrast, the photon counting method has a very 
elegant generalization to the multimode case: after applying the displacement to each of the in- 
volved modes, the Wigner function at the selected point is simply given by the average parity of 
the total number of detected photons. Moreover, the dichotomic outcome of such a measurement 
provides a novel way of testing quantum nonlocality exhibited by correlated states of optical 
radiation 



Chapter 7 

Statistical uncertainty in photodetection 
measurements 



In this chapter, we shall study the problem of statistical uncertainty from a more general point of 
view. Over recent years, the set of tools for measuring quantum statistical properties of optical 
radiation has substantially enlarged. In addition to the techniques discussed in this thesis: dou- 
ble homodyne detection, optical homodyne tomography, and the direct method for measuring 
quasidistribution functions, other novel schemes have been proposed p^, B51 EB, E7|, FHp. The 



quantum optical "toolbox" for measuring light contains now experimentally established schemes 
for reconstructing various representations of its quantum state: the Q function, the Wigner func- 
tion, and the density matrix in the quadrature and the Fock bases. A device that is used in most 
of quantum optical schemes to convert the quantum signal to a macroscopic level is the photode- 
tector. Thus, photodetection is a basic ingredient of quantum optical measurements. 

The quantum state can be characterized using various representations: quasidistribution func- 
tions or a density matrix in a specific basis. From a theoretical point of view, all these forms are 
equivalent. Each representation contains complete information on the quantum state, and they 
can be transformed from one to another. Any observable related to the measured systems can be 



evaluated from an arbitrary representation using an appropriate expression [pOQ. 

This simple picture becomes much more complicated when we deal with real experimental 
data rather than analytical formulae. Each quantity determined from a finite number of exper- 
imental runs is affected by a statistical error. Consequently, the density matrix or the quasidis- 
tribution function reconstructed from the experimental data is known only with some statistical 
uncertainty. This uncertainty is important when we further use the reconstructed information to 
calculate other observables or to pass to another representation. The crucial question is, whether 
determination of a certain representation with sufficient accuracy guarantees that arbitrary ob- 
servable can be calculated from these data with a reasonably low statistical error. If this is not 
the case, the reconstructed information on the quantum state of the measured system turns out 
to be somewhat incomplete. Furthermore, transformation between various representations of the 
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quantum state becomes a delicate matter. 

These and related problems call for a rigorous statistical analysis of quantum optical mea- 
surements. We shall provide here complete statistical description of a measurement of quantum 
observables in optical schemes based on photodetection. This approach fully characterizes sta- 
tistical properties of quantities determined in a realistic measurement from a finite number of 
experimental runs. It can be applied either to determination of a single quantum observable, 
or to the reconstruction of the quantum state in a specific representation. Within the presented 
framework we study, using a simple example, the completeness of the experimentally recon- 
structed information on the quantum state. We demonstrate the pathological behavior suggested 
above, when the reconstructed data cannot be used to calculate certain observables due to rapidly 
exploding statistical errors. Our example is motivated by the direct scheme for measuring qua- 
sidistribution functions of light. 

Compared to previous works on statistical noise, our attention will be focused here on two im- 
portant statistical aspects of quantum state measurement. First, we shall go beyond the minimum 



second-order treatment of the statistical error [ pO| , pTj , p2| ] and provide, in a closed mathematical 
form, a complete statistical description of quantum observables determined from realistic mea- 
surements. This result, derived from the first principles without using any approximations, is 
an exact analytical solution to the problem which so far has been approached only by means of 
Monte Carlo simulations [HR B^. The second problem discussed here will be the feasibility of 



reconstructing the quantum state from realistic, finite data collected in a specific experimental 
scheme. It is generally believed that a sufficient condition for successful reconstruction is the 
existence of the co variance matrix with finite elements [|15], |16[ [T^]. We point out that this belief 
misses an important issue resulting from statistical uncertainty. Suppose we have reconstructed 
a family of observables with finite statistical variances. In order to be sure that the reconstructed 
information on the quantum state is accurate and complete, we should ask the following question: 
can we always use these observables to evaluate any quantum property of the measured system 
that can be expressed in terms of the reconstructed family? We shall give a clear negative answer 
to this question, based on a detailed discussion of carefully selected counterexamples. Although 
all observables involved in these examples are represented by bounded, well-defined operators, 
statistical fluctuations are shown to be arbitrarily huge, when the reconstructed family is further 
used to evaluate certain quantum expectation values. This singular behaviour is quantitatively ex- 
plained as a result of strong statistical correlations between observables reconstructed from the 
same sample of experimental data. Consequently, the discussed examples clearly demonstrate 
that statistical properties of reconstructed observables depend in an essential way on a specific 
experimental scheme, which effectively limits available information on the measured system. 

This chapter is organized as follows. The starting point of our analysis is the probability 
distribution of obtaining a specific histogram from N runs of the experimental setup. This basic 
quantity determines all statistical properties of quantum observables reconstructed from a finite 
sample of experimental data. We characterize these properties using the generating function. 
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for which we derive an exact expression directly from the probability distribution of the exper- 
imental outcomes. These general results are presented in Sec. [77T] . Then, in Sec. [7!^ we use 
the developed formalism to discuss the reconstruction of the photon number distribution of a 
single light mode, and its subsequent utilization to evaluate the parity operator 11. We consider 
two experimental schemes: direct photon counting using an imperfect detector, and homodyne 
detection with random phase. In both the cases we find that the evaluation of the parity operator 
from the reconstructed photon statistics is a very delicate matter. For photon counting of a ther- 
mal state, we show that neither the statistical mean value of 11 nor its variance have to exist when 
we take into account arbitrarily high count numbers. For random phase homodyne detection, the 
statistical error of the parity operator is an interplay of the number of runs A^ and the specific 
regularization method used for its evaluation. The example of the parity operator illustrates diffi- 
culties related to the transformations between various experimentally determined representations 
of the quantum state, as the parity operator yields, up to a multiplicative constant, the Wigner 



function at the phase space origin. Finally, in Sec. ^3 ^^ discuss consequences of statistical 
uncertainty in the measurements of the quantum state. 

7.1 Statistical analysis of experiment 

In photodetection measurements, the raw quantity delivered by a single experimental run is the 
number of photoelectrons ejected from the active material of the detectors. The data recorded 
for further processing depends on a specific scheme. It may be just the number of counts on 
a single detector, or a difference of photocounts on a pair of photodetectors, which is the case 
of balanced homodyne detection. It may also be a finite sequence of integer numbers, e.g. for 
double homodyne detection. We will denote in general this data by n, keeping in mind all the 
possibilities. 

The experimental scheme may have some external parameters 9, for example the phase of 
the local oscillator in homodyne detection. The series of measurements are repeated for various 
settings 9i of these parameters. Thus, what is eventually obtained from the experiment, is a 
set of histograms {kn}ei, telling in how many runs with the settings 9i the outcome n has been 
recorded. We assume that for each setting the same total number of A^ runs has been performed. 

The theoretical probability Pn{9i) of obtaining the outcome n in a run with settings 9i is 
given by the expectation value of a positive operator- valued measure (POVM) Pn{9i) acting in 
the Hilbert space of the measured system. The reconstruction of an observable A is possible, if it 
can be represented as a linear combination of the POVMs for the settings used in the experiment: 

i n 

where an{9i) are the kernel functions. The above formula allows one to compute the quantum 
expectation value {A) from the probability distributions ]9„(^j) = (p„(^i)). 
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This theoretical relation has to be applied now to the experimental data. The simplest and 
the most commonly used strategy is to estimate the probability distributions Pn{Oi) by experi- 
mental relative frequencies {kn/N)0.. The relative frequencies integrated with the appropriate 
kernel functions yield an estimate for the expectation value of the operator A. For simplicity, we 
will denote this estimate just by A. Thus, the recipe for reconstructing the observable A from 
experimental data is given by the counterpart of Eq. (|77T|): 

in. ^ '^ "i 



We will now analyse statistical properties of the observable A evaluated according to Eq. ( [7. 2D 
from data collected in a finite number of experimental runs. Our goal is to characterize the 
statistical distribution w{A) defining the probability that the experiment yields a specific result 
A. The fundamental object in this analysis is the probability V{{kn};9) of obtaining a specific 
histogram {kn} for the settings 9. In order to avoid convergence problems, we will restrict the 
possible values of n to a finite set by introducing a cut-off. The probability V{{kn}; 6) is then 
given by the multinomial distribution [ |5D| ] : 

p({^.};^)= ,^_^,,_), [i-E>)) n'j^w«)i'". (7.3) 

where prim in sums and products denotes the cut-off. This distribution describing experimental 
histograms is derived from an assumption that all the detection events are statistically indepen- 
dent. 

Let us first consider a contribution Ai to the observable A calculated from the histogram Oi'. 

A. = Yl^n{e,){^) ■ (7-4) 

Its Statistical distribution w{Ai] Oj) is given by the following sum over all possible histograms 
that can be obtained from N experimental runs: 

w{Af,ei) = J2v{{kr,};e,)6 I A, - lj^an{e.)kA . (V.S) 

jfenl \ n / 



Equivalently, the statistics of Ai can be characterized by the generating function w{X; Oi) for the 
moments, which is the Fourier transform of the distribution w{Ai] 6i): 



w{X;e,) = / dAe'^'^'wiA^Oi) 



J2 '^a^nh e,) exp ( !A^'a„(^,)A;„ ) . (7.6) 

{kn} 
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An easy calculation yields the explicit form of the generating function: 



w{X;d, 



l + Y^Pn{e^){e 



N 



i\a„{ei)/N 



(7.7) 



The observable A is obtained via summation of the components Ai corresponding to all settings 
of the external parameters 6i. As these components are determined from disjoint subsets of the 
experimental data, they are statistically independent. Consequently, the generating function w{\) 
for the moments of the observable A is given by the product: 



w 



(A) = UAe'^MA) = \{w{\:ei) 

i \ n / 



N 



(7.8) 



This expression contains the complete statistical information on determination of the observable 
A from a finite number of runs of a specific experimental setup. The measuring apparatus is 
included in this expression in the form of a family of POVMs p„(6'j). The quantum expectation 
value of these POVMs over the state of the measured system yields the probability distributions 
Pn{Oi)- Finally, the coefficients an{9i) are given by the computational recipe for reconstructing 
the observable A from the measured distributions. The analytical expression for the generating 
function given in Eq. (|7.8|), derived from the exact description of raw experimental outcomes, 
provides a complete characterization of statistical fluctuations in realistic measurements of quan- 
tum observables. Let us note that in general the generating function w{\) cannot be expressed by 
the operator A alone. Both the POVMs Pn{Oi) and the kernel functions an{6i) enter Eq. ( p.8[ ) in 
a nontrivial way, which clearly shows that statistical properties of the reconstructed observable 
depend essentially on the specific measurement scheme. 

The basic characteristics of statistical properties of the observable A is provided by the mean 
value E(A) and the variance Var(y4). These two quantities can be easily found by differentiating 
the logarithm of the generating operator: 



E{A) :- 
Var(A) := 



-— logw(A) 
t dA 



A=0 



y^y^Qn(6'»)pn 



(7.9) 



A=0 



1 

N 



Y,Y^<^nmVn{e{) - J2 [Y!a^{0^)Pn{ 



(7.10) 



The statistical error is scaled with the inverse of the square root of the number of runs A^. Let us 
note that the second component in the derived formula for Var(A) differs from that used in the 
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discussions of homodyne tomography in Refs. P^ , |5^ ], where it was equal just to [E(yl)]^. This 
difference results from different assumptions about the local oscillator phase: in Refs. [ ^ |5^ 
it was considered to be a uniformly distributed stochastic variable in order to avoid systematic 
errors, whereas we have assumed that the number of runs is fixed for each selected setting of the 
external parameters. 

The goal of quantum state measurements is to retrieve the maximum amount of informa- 
tion on the quantum state available from the experimental data. Therefore the experimental 
histograms are usually processed many times in order to reconstruct a family of observables 
characterizing the quantum state. Of course, quantities determined from the same set of exper- 
imental data are not statistically independent, but in general exhibit correlations. The analysis 
presented above can be easily extended to the evaluation of any number of observables from the 
same sample of experimental data. If we restrict our attention to the basic, second-order char- 
acterization of these correlations, it is sufficient to discuss simultaneous determination of two 
observables. Let us suppose that in addition to A, another observable B has been calculated 
from the histograms {knjei according to the formula: 



B 



EE''"( 



N 



(7.11) 



The generating function w{\, jj) corresponding to the joint probability distribution w{A, B) can 
be found analogously to the calculations presented above. The final result is: 



w(A,/i) = fdAdBe'^'^+'''^w{A,B) 

i \ n 



N 



iXa„{9,)/N+ifj.bn{9^)/N _ -^\ 



(7.12) 



The covariance between the experimentally determined values of A and B is given by: 



Cov{A,B) 



1 d^ 
P dAd/i 



logw{\,n) 



X,fj.=0 



^E Y^cin{ei)bn{e.)Pn{e^) - (Y!an{e.)Pn{ei)] (^'6„.(^.)p, 




(7.13) 



The covariance can be normalized to the interval [—1,1] using Var(A) and Var(i?), which yields 
the correlation coefficient for the pair of observables A and B: 



Corr{A,B) :-- 



Co\{A, B) 
v/Var(A)Var(5) ' 



(7.14) 
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This quantity defines whether the statistical deviations of A and B tend to have the same or 
opposite sign, which corresponds respectively to the positive or negative value of Corr(y4, B). 

We have assumed that the histograms A;„ have been measured for a finite number of external 
parameters settings 9i, which is always the case in an experiment. However, in some schemes the 
measurement of histograms is in principle necessary for all values of a continuous parameter. For 
example, in optical homodyne tomography the full information on the quantum state is contained 
in a family of quadrature distributions for all local oscillator phases. Restriction to a finite set of 
phases introduces a systematic error to the measurement [ |5D| , [7D| ]. 

7.2 Phase-insensitive detection of a light mode 

We will now apply the general formalism developed in the preceding section to the reconstruc- 
tion of phase-independent properties of a single light mode. The basic advantage of this exem- 
plary system is that it will allow us to discuss, in a very transparent way, pathologies resulting 
from the statistical uncertainty. All phase-independent properties of a single light mode are 
fully characterized by its photon number distribution p,^. Therefore, it is sufficient to apply a 
phase-insensitive technique to measure the photon statistics of the field. We will consider two 
measurement schemes that can be used for this purpose: direct photon counting and random 
phase homodyne detection. 

The photon number distribution is given by the expectation value of a family of projection 
operators p^ = |z/) (z/|, where |z/) is the z/th Fock state. In principle, knowledge of this distribution 
enables us to evaluate any phase-independent observable related to the measured field. A simple 
yet nontrivial observable, which we will use to point out difficulties with the completeness of the 
reconstructed information on the quantum state, is the parity operator: 

n = J](-irk)(H- (V.15) 

iy=0 

This operator is bounded, and well defined on the complete Hilbert space of a single light mode. 
Its expectation value is given by the alternating series of the photon number distribution: 

oo 

(n) = j2{-inp.), (7.16) 

iy=0 

which is absolutely convergent for any quantum state. Therefore, any pathologies connected to 
its determination from experimental data, if there appear any, cannot be ascribed to its singular 
analytical properties. 
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7.2.1 Direct photon counting 

First, we will consider the reconstruction of phase insensitive properties of a single light mode 
from data measured using a realistic, imperfect photodetector. The positive operator- valued mea- 
sure Pn describing the probability of ejecting n photoelectrons from the detector is given by [[71]]: 

Pn = ■ j — exp(-?7a'a) :, (7.17) 

where a is the annihilation operator of the light mode, and r] is the quantum efficiency of the 
photodetector. In the limit r^ — ;► 1 we get directly p„ = \n){n\. In a general case, the probability 
distribution for the photoelectron number is related to the photon statistics via the Bernoulli 
transformation. This relation can be analytically inverted [JTSl], which yields the expression: 

oo 



n=0 



where the kernel functions ruJ are given by: 



0, n < u, 



The inversion formula has a remarkable property that p^, depends only on the "tail" of the pho- 
tocount statistics for n > u. It has been shown that the inverse transformation can be applied to 
experimentally determined photocount statistics for an arbitrary state of the field, provided that 



the detection efficiency is higher that 50% [[ISp. 



We will now discuss statistical properties of the photon number distribution determined by 



photon counting within the general framework developed in Sec. ^TT]- The case of perfect de- 
tection is trivial for statistical analysis. Therefore we will consider nonunit detection efficiency, 
which is numerically compensated in the reconstruction process using the inverse Bernoulli trans- 
formation according to Eq. ( [7.18p . For all examples presented here, the efficiency is 77 = 80%, 
which is well above the 50% stability limit. 



In Fig. ^J] we depict the reconstructed photon number distributions for a coherent state, a 
thermal state, and a squeezed vacuum state. The mean values E(p,^) along with their statistical 
errors [Var(p,^)]^/^ are compared with Monte Carlo realizations of a photon counting experiment, 
with the number of runs A^ = 4000. It is seen that for a thermal state and a squeezed vacuum 
state, the statistical error of the probabilities pi^ grows unlimitedly with the photon number u. 
However, any experimental histogram obtained from a finite number of runs ends up for a cer- 
tain count number, and therefore the reconstructed photon statistics is zero above this number. 
An important feature that is evidently seen in the Monte Carlo simulations, are correlations be- 
tween the consecutive matrix elements. The reconstructed photon number distribution clearly 
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Figure 7.1: Reconstruction of the photon number distribution from photon counting for (a) a coherent 
state with (n) = 4, (b) a squeezed vacuum state with (n) = 1, and (c) a thermal state with (n) = 2, from 
N = 4000 runs in each case. Monte Carlo simulations of a photon counting experiment, depicted with 
points, are compared with exact values (solid lines), with the statistical errors [Var(p,y)]^/^ marked as grey 
areas. The detection efficiency is 77 = 



exhibits oscillations around the true values. This property can be quantified using the correlation 



coefficient defined in Eq. ( [7.14| ), which we plot for all three states in Fig. \1 .% For large z/'s, 
Corr(p,^,P;^+i) is close to its minimum allowed value —1, which acknowledges that statistical 
correlations are indeed significant. 

These correlations affect any quantity computed from the reconstructed photon number dis- 
tribution. The parity operator is here a good example: since in Eq. ( |7.16D we sum up consecutive 
Pj,'s with opposite signs, their statistical deviations do not add randomly, but rather contribute 
with the same sign. Consequently, the statistical error of the evaluated parity operator may be 
huge. It is therefore interesting to study this case in detail. In principle, we could obtain statisti- 
cal properties of the reconstructed parity using the covariance matrix for the photon distribution. 
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Figure 7.2: The correlation coefficient between the consecutive density matrix elements Corr:{py , p^^i] 
depicted for the coherent state (o), the squeezed state (D), and the thermal state (*) from Fig. Ol. 



However, it will be more instructive to express the parity directly in terms of the photocount 
statistics, and then to apply the statistical analysis to this reconstruction recipe. This route is 
completely equivalent to studying evaluation of the parity via the photon number distribution, as 
all transformations of the experimental data, which we consider here, are linear. 
A simple calculation combining Eqs. ( |7.15| ) and ( [7.18[ ) shows that: 

n = 5Z 1-- Pn, (7.20) 



n=0 



where we have introduced in the upper summation limit a cut-off parameter K for the photocount 
number. This formula clearly demonstrates pathologies related to the determination of the parity 
operator. For any rj < 1, the factor (1 — 2/77)" is not bounded, which makes the convergence of 
the whole series questionable in the limit i^ — >^ 00. Of course, for an experimental histogram the 
summation is always finite, but the exploding factor amplifies contribution from the "tail" of the 
histogram, where usually only few events are recorded, and consequently statistical errors are 
significant. 

Let us study these pathologies more closely using examples of a coherent state and a thermal 
state. For a coherent state |a), both the expressions for E(n) and Var(n) are convergent with 
i^ ^ 00. However, the variance, given by the formula: 



Var(n' 



coh\ 



1 



exp 



4(1 - r/) 

7] 



\a\ 



exp(— 4|a| 



(7.21) 



grows very rapidly with the coherent state amplitude a, when the number of runs N is fixed. 
For a thermal state with the average photon number n, the matter becomes more delicate. When 
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Figure 7.3: Determination of the parity operator for (a) coherent states and (b) thermal states with the in- 
creasing average photon number, assuming the photodetector efficiency t] = 80%. Each square represents 
the parity evaluated from Monte Carlo simulated photon statistics with N = 4000 runs for a given average 
photon number. The solid lines and the grey areas depict the mean value E(n) and the error [Var(n)]^/^. 
For thermal states, the variance diverges to infinity when n > 1, and the mean value E(n) does not exist 
above n > 2.5. 



K ^ oo, the series ( [7 .201 ) is convergent only for n < 1/[2(1 — r/)], which for i] = 80% gives just 
2.5 photons. Even when the mean value exists, the variance is finite only for n < r]/[4:{l — t])] 
and equals: 



1 



Var(n''') = - 



V 



(7.22) 



We illustrate these results with Fig. [73| , depicting Monte Carlo simulations for various average 
photon numbers. For coherent states, statistical fluctuations can in principle be suppressed by 
increasing the number of runs. For thermal states, the situation is worse: when n > I, the 
variance cannot even be used as a measure of statistical uncertainty. 

Let us recall that the bound r] > 50% for the stability of the inverse Bernoulli transformation 
is independent of the state to be measured. It has been obtained from the requirement that in 
the limit i^ ^ oo both ^{pu) and Var(p,^) should converge [[TSl]. The example with the parity 
operator clearly shows, that the condition i] > 50% does not guarantee that the reconstructed 
photon number distribution can be safely used to determine an arbitrary well-behaved phase 
independent observable. Thus, imperfect detection is inevitably connected with some loss of the 
information on the measured quantum state. 

We have noted that as long as finite, experimental data are concerned, evaluation of ob- 
servables via intermediate quantities is equivalent to expressing them directly in terms of the 
measured probability distributions. One might try to circumvent the 7] > 50% bound for recon- 
structing the photon statistics by applying the inverse Bernoulli transformation in two or more 
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steps, and compensating in each step only a fraction of its inefficiency. Of course, such a strat- 
egy must fail, as for any finite sample of experimental data this treatment is equivalent to a single 
transformation which is unstable. In many-step processing this instability would be reflected in 
increasing correlations and statistical errors exploding to infinity. 

7.2.2 Random phase homodyne detection 



As we discussed in Sec. [3^ , random phase homodyne detection is a recently developed technique 
for measuring phase-independent properties of optical radiation, which goes beyond certain lim- 
itations of plain photon counting p5|]. Data recorded in this scheme is the difference of counts 



on two photodetectors measuring superposition of the signal field with a strong coherent local 
oscillator. The count difference is rescaled by the local oscillator amplitude, and the resulting 
stochastic variable x can be treated with a good approximation as a continuous one. The pho- 
ton number distribution is reconstructed from the random phase homodyne statistics pni^) by 
integrating it with pattern functions fuix) [ p6| , p7| , ^: 

/oo 
dxf^{x)pTi{x). (7.23) 

A convenient method for numerical evaluation of the pattern functions has been described in 
Ref. [0]. 

Let us now discuss statistical properties of the homodyne scheme in its discretized version 
used in experiments, when the rescaled count difference is divided into finite width bins. As the 
local oscillator phase is random, the setup has no controllable parameters, and the statistics of 
the observables is fully determined by P7^(a;). Statistical errors of the density matrix in the Fock 
basis reconstructed via homodyne detection have been studied in Ref. [^. Here we will focus 
our attention on statistical correlations exhibited by the diagonal density matrix elements, and 
their further utilization for evaluating phase-independent observables. 

We will consider the unit detection efficiency i] = I, with no compensation in the processing 
of the experimental data. This is the most regular case from the numerical point of view. When 
1] < I and the compensation is employed, the statistical errors are known to increase dramatically 



[pOQ. In Fig. ^4| we depict the homodyne reconstruction of the photon number distribution for 
the three states discussed in the previous subsection. For large u, the statistical errors tend to 
a fixed value >/2jN, which has been explained by D'Ariano et al. using the asymptotic form 



of the pattern functions [Q. Again, Monte Carlo simulations suggest that the reconstructed 
density matrix elements are correlated, which is confirmed by the correlation coefficient for the 
consecutive photon number probabilities, plotted in Fig. [73| . A simple analytical calculation 
involving the asymptotic form of the pattern functions shows, that for large u this coefficient 
tends to its minimum value —1. 

One may now expect that no subtleties can be hidden in using the reconstructed photon 
number distribution to evaluate the parity operator according to Eq. ( [7.16D . However, let us 
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Figure 7.4: Random phase homodyne reconstruction of the photon number distribution for (a) the coherent 



state (b) the squeezed state and (c) the thermal state, all states with the same photon numbers as in Fig. \Ll 
The range of the homodyne variable is restricted to the interval — 6 < x < 6 divided into 1200 bins. The 
simulated homodyne statistics is obtained from A^ = 4 • 10^ Monte Carlo events. 



recall that the parity operator is equal, up to a multiplicative constant, to the Wigner function 
at the phase space origin. The Wigner function is related to the homodyne statistics via the 
inverse Radon transformation, which is singular. In particular, applying this transformation for 
the phase space point (0, 0) we obtain the following expression for the parity operator in terms 
of the homodyne statistics [ [7^ : 



(n) 



dxpTzix)-—P-, 
ax X 



(7.24) 



where P denotes the principal value. Due to the singularity of the inverse Radon transform, its 
application to experimental data has to be preceded by a special filtering procedure. This feature 
must somehow show up, when we evaluate the parity operator from the reconstructed photon 
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Figure 7.5: The correlation coefficient for the consecutive photon number probabilities for the coherent 
state (o), the squeezed state (D) and the thermal state (*). 



statistics. In order to analyse this problem in detail let us discuss evaluation of the truncated 
parity operator fix from a finite part of the photon number distribution: 



(n^) = 5^(-l)^(p.) = / gK{x)p{x), 

n J —CXI 



where 



u=0 



9k{x) 



K 

E 

v=0 



i-^rux) 



(7.25) 



(7.26) 



can be considered to be a regularized kernel function for the parity operator. In Fig. ^^ we plot 
this function for increasing values of the cut-off parameter K. It is seen that the singularity of 
the kernel function in the limit fC ^ oo is reflected by an oscillatory behaviour around x = 
with growing both the amplitude and the frequency. This amplifies the statistical uncertainty of 
the experimental homodyne data. In Fig. ^^ we show determination of the parity operator for 
the three states discussed before, using increasing values of the cut-off parameter K. Though we 
are in the region where the true photon number distribution is negligibly small, addition of sub- 
sequent matrix elements increases the statistical error in an approximately linear manner. This 
is easily understood, if we look again at the reconstructed photon number distributions: increas- 
ing K by one means a contribution of the order of ^/2/N added to the statistical uncertainty, 
and, moreover, these contributions tend to have the same sign due to correlations between the 
consecutive matrix elements. Thus, determination of the parity operator from homodyne statis- 
tics requires an application of a certain regularization procedure. It may be either the filtering 
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Figure 7.6: Regularized kernel functions for the parity operator gxi^) evaluated as a finite sum of the 
Fock states pattern functions, for increasing values of the cut-off parameter K. 



used in tomographic back-projection algorithms, or the cut-off of the photon number distribution. 
The statistical uncertainty of the final outcome is eventually a result of an interplay between the 
number of experimental runs and the applied regularization scheme. 

Finally, let us briefly comment on the compensation for the nonunit efficiency of the homo- 
dyne detector. First, one might think of applying a two-mode inverse Bernoulli transformation 
directly to the joint count statistics on the detectors. However, it is impossible in the homodyne 
scheme to resolve contributions from single absorbed photons due to high intensity of the de- 
tected fields. The inverse Bernoulli transformation has no continuous limit, as consecutive count 
probabilities are added with opposite signs. Nevertheless, the nonunit detection efficiency can 
be taken into account in the pattern functions p7| , ^Bp . In this case the statistical errors increase 
dramatically, and explode with z/ — > oo, which makes determination of the parity operator even 
more problematic. This is easily understood within the phase space picture: the distributions 
measured by an imperfect homodyne detector are smeared-out by a convolution with a Gaus- 
sian function [ p5| , p3| ] . Evaluation of the parity operator, or equivalently, the Wigner function 
at the phase space origin requires application of a deconvolution procedure, which enormously 
amplifies the statistical error [[73]] . 



7.3 Consequences for quantum state measurement 



We have presented a complete statistical analysis of determining quantum observables in optical 
measurement schemes based on photodetection. We have derived an exact expression for the 
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Figure 7.7: Reconstruction of the truncated parity operator IIk for the coherent state (o), the squeezed 
state (D), and the thermal state (*) with various values of the cut-off parameter K, using the same Monte 
Carlo homodyne statistics as in Fig. [7.4[ The simulations are compared with the corresponding mean 
values E{IIk) and errors [Var(n/^)]^'^ plotted as solid lines surrounded by grey areas. 



generating function characterizing statistical moments of the reconstructed observables, which, 
in particular, provides formulae for statistical errors and correlations between the determined 
quantities. These general results have been applied to the detection of phase-independent prop- 
erties of a single light mode using two schemes: direct photon counting, and random phase 
homodyne detection. This study has revealed difficulties related to the completeness of the re- 
constructed information on the quantum state: in some cases the parity observable, which is a 
well-behaved bounded operator, effectively cannot be evaluated from the reconstructed data due 
to the exploding statistical error. 

We have recalled that the parity operator is directly related to the value of the Wigner function 
at the phase space origin. Thus, our example can also be interpreted as a particular case of the 
transformation between two representations of the quantum state: in fact, we have considered 
evaluation of the Wigner function at a specific point (0, 0) from the relevant elements of the 
density matrix in the Fock basis. Therefore, our discussion exemplifies subtleties related to the 
transition between various quantum state representations, when we deal with data reconstructed 
in experiments. Though a certain representation can be determined with the statistical uncertainty 
which seems to be reasonably small, it effectively cannot be converted to another one due to 
accumulating statistical errors. 

The presented study suggests that the notion of completeness in quantum state measurements 
should inherently take into account the statistical uncertainty. From a theoretical point of view, 
the quantum state can be characterized in many different ways which are equivalent as long as 
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expectation values of quantum operators are known with perfect accuracy. In a real experiment, 
however, we always have to keep in mind the specific experimental scheme used to perform 
the measurement. This scheme defines statistical properties of the reconstructed quantities, and 
may effectively limit the available information on the quantum state. Determination of a family 
of observables does not automatically guarantee the feasibility of reconstructing the expectation 
value of an arbitrary well behaved operator. Reconstruction of any observable should be preceded 
by an analysis how significantly the final result is affected by statistical noise corrupting the raw 
experimental data. 
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Chapter 8 
Conclusions 



In this thesis, we have developed a direct method for measuring quasidistribution functions of 
light. Quasidistribution functions provide a complete characterization of the quantum state in the 
form which is analogous to a classical phase space distribution. We have shown that quasidis- 
tributions, in particular the Wigner function, can be determined using an optical scheme based 
on photon counting. We have reported an experimental realisation of this scheme, and presented 
measurements of the Wigner function for several classical-like states. 

There is a wide range of problems arising from a closer look at various aspects of the proposed 
scheme. We have discussed the deleterious effects of imperfect detection, and shown that they 
cannot be in general compensated in numerical processing of the experimental data. Further, 
we have developed a general multimode theory of the scheme, and extended the principle of 
the measurement to multimode quasidistributions. The developed multimode theory has been a 
useful tool in studying the role of experimental imperfections. We have also seen that in order 
to analyse the feasibility of the measurement, it is necessary to take into account the statistical 
uncertainty. We have derived estimates for the statistical error, and shown that in some cases the 
statistical noise effectively limits available information on the quantum state. 

The field of quantum state measurement is continuously developing. Currently, an interesting 
direction of research is related to fundamental aspects of retrieving information on the quantum 
state. Generally, our source of information on the quantum state is the result of a measurement 
performed on an ensemble of a finite number of copies. From these data we want to infer the 
quantum state of the ensemble. The inference can be performed in different ways, based on 
various statistical methodologies. Throughout this thesis, we have used the linear approach, 
where quantum probability distributions are estimated by relative frequency histograms. This is 
the most straightforward and so far the most commonly used strategy, but it is not necessarily 
the most efficient one. Recently, other approaches to data processing have been proposed, based 
on the maximum entropy principle [[7^], the least-squares inversion [|75|], and the maximum- 
likelihood estimation [ [7^ , |57| , [77| , [7S| ]. These methods enhance the amount of information which 
can be obtained from a realistic measurement. From a more fundamental point of view, one 
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may ask what is the ultimate bound for estimating the quantum state from finite ensembles. This 
question can be reformulated as a problem of designing the optimal quantum state measurement 
on finite ensembles. Such a problem is in general extremely difficult, and so far it has been 
studied only in selected yet highly nontrivial cases [^, ^ |81]]. The optimal strategy has been 
shown to be of the form of a collective measurement performed jointly on all the available copies. 
The research on quantum state measurement significantly contributes to the understanding of 
the foundations of quantum theory. Quantum mechanics still unveils surprising consequences of 
its principles. An important recent example is the observation that quantum mechanics opens up 
completely new ways of processing, storing, and transmitting information. Exploration of these 
these possibilities can result in novel information technologies, such us quantum computing [^ 
and quantum cryptography [jS^, |S^]. Some issues in quantum state measurement are shared with 
the quantum information theory, for example methods of quantum estimation and the effect of 
decoherence. On the other hand, quantum information technologies need carefully prepared 
quantum systems, and measuring the quantum state is here an indispensable diagnostic tool. 
With no doubt, the field of quantum state measurement plays a prominent role at the forefront of 
contemporary physics. 
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(a) Vacuum state 



(b) Coherent state 



(c) Phase diffused coherent state 







The measured Wigner functions for (a) the vacuum, (b) a weak coherent state with approximately one photon, and (c) a phase diffused coherent state. The photon 
statistics was collected on a polar grid spanned by 20 amplitudes, and 50 phases for the plots (a) and (b), or 40 phases for the plot (c). The duration of a single 
counting interval was 40 fis for (a) and (b) and 30 ^s for (c). The measurements were performed for slightly different laser intensities, and the radial coordinate for 
each of the graphs was scaled separately. 
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